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1.  INTRODUCTION 


A  metamaterial  is  a  material  that  gains  its  properties  from  its  structure  rather  than 
from  its  constituent  material  phases.  Metamaterials  with  negative  refractive  index  is  one 
of  the  most  active  topics  among  the  research  community.  The  idea  of  it  traces  back  to 
1968  when  Veselago  (1968)  postulated  a  theory  for  possible  materials  having  negative 
electric  permittivity  (s)  and  magnetic  permeability  (p),  and  hence  resulting  in  a  negative 
refractive  index.  Three  decades  later,  Smith  et  al.  (2000)  proposed  a  design  of  making  an 
artificial  material  that  exhibits  a  frequency  region  in  the  microwave  regime  with 
simultaneously  negative  permittivity  and  negative  permeability.  They  demonstrated  both 
by  numerical  simulation  and  experiment  the  unusual  physical  phenomena.  In  the  same 
year,  Pendry  (2000)  proposed  a  practical  application  of  this  type  of  electromagnetic 
metamaterials:  a  perfect  lens,  capable  of  imaging  objects  and  fine  structures  that  are 
much  smaller  than  the  wavelength  of  light.  As  a  result  of  this  development  of  novel 
electromagnetic  metamaterials,  researchers  extended  the  same  concept  to  the  acoustic 
counterparts  immediately.  For  the  negative-index  acoustic  metamaterials,  basically  it  is  to 
find  materials  possessing  negative  mass  density  and  negative  modulus  because  these  two 
material  properties  govern  the  wave  motions  in  acoustic  (or  elastic)  materials.  However, 
in  contrast  to  the  natural  negative  electric  permittivity,  no  natural  material  has  been  found 
to  possess  negative  mass  density  or  negative  modulus.  This  difficulty  brings  about  some 
barriers  to  the  realization  of  acoustic  counterparts. 

Some  recent  developments  have  theoretically  shown  that  acoustic/elastic 
metamaterials  could  be  realized.  Liu  et  al.  (2000)  investigated  a  new  type  of  composite 
fabricated  using  lead  balls  coated  with  silicon  rubber  in  the  epoxy  matrix.  Around  the 
local  resonance  frequencies,  the  composite  behaves  as  a  material  with  effective 
negativity.  Such  an  anomalous  response  of  the  system  is  the  consequence  of  the  coupling 
of  the  long-wavelength  elastic  wave  in  the  host  medium  with  the  localized  vibration 
motion  of  the  spheres.  Recently,  Milton  and  Willis  (2007)  and  Huang  et  al.  (2009)  have 
presented  a  theoretical  foundation  for  the  concept  of  negative  effective  mass. 

The  field  of  metamaterials  is  still  in  its  infancy.  Further  studies  are  needed  to  fully 
understand  the  complex  behavior  of  these  materials.  More  designs  of  microstructures  are 
desired  to  produce  new  metamaterials  with  new  properties  and  some  of  the  unusual 
behavior  of  metamaterials  predicted  by  theoretical  models  need  to  be  verified 
experimentally.  More  importantly,  practical  applications  of  metamaterials  should  be 
explored. 


2.  OBJECTIVE 

The  purpose  of  this  research  is  to  design,  fabricate,  and  study  new  acoustic 
metamaterials  with  a  negative  effective  mass  density  and  modulus.  Through  this  study,  it 
is  anticipated  that  the  implications  to  material  and  structural  applications  of  the  negative 
effective  mass  and  modulus  of  metamaterials  can  be  fully  understood.  Moreover,  more 
microstructures  for  metamaterials  will  be  designed  and  stidied.  The  3D  printing 
technology  will  be  utilized  to  produce  these  complicated  microstructures  and 
metamaterial  specimens  for  experiments. 
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3.  APPROACH 


The  key  feature  of  a  metamaterial  that  may  exhibit  a  negative  effective  mass  density 
and  modulus  takes  the  form  of  a  continuum  containing  local  resonators.  This 
configuration  yields  a  local  resonance  frequency,  a  feature  essential  to  the  desired 
metamaterial.  In  view  of  this,  we  focus  the  effort  on  the  design  of  new  resonator 
microstructures,  perform  modeling  and  simulations  of  the  dynamic  behavior  of  these 
materials,  conduct  experiments  to  verify  the  model  predictions,  and  identify  applications 
of  these  metamaterials. 


4.  Acoustic  Metamaterial  with  Extreme  Young’s  Moduli 

4.1  Introduction 

The  subject  of  metamaterials  has  recently  attracted  much  attention  of  researchers. 
Generally  speaking,  metamaterials  are  considered  a  new  class  of  composites  with  man¬ 
made  subwavelength  microstructures  that  exhibit  unusual  behavior  not  readily  observed 
in  natural  materials.  Veselago  (1968)  first  postulated  a  theory  for  possible  materials 
having  negative  electric  permittivity  and  magnetic  permeability,  and  hence  resulting  in  a 
backward  wave  and  a  negative  refractive  index.  This  concept,  however,  did  not  stimulate 
much  interest  among  researchers  at  that  time.  Three  decades  later.  Smith  et  al.  (2000) 
investigated  the  possibility  of  making  left-handed  electromagnetic  metamaterials.  Pendry 
(2000)  proposed  a  perfect  lens  made  of  left-handed  metamaterial  that  would  lead  to 
practical  applications.  Shelby  et  al.  (2001)  demonstrated  the  experimental  verification  of 
a  negative  index  of  refraction.  Since  then,  numerous  works  regarding  metamaterials  have 
been  engaged  and  great  numbers  of  publications  are  poured  in. 

Motivated  by  the  mathematical  analogy  between  acoustic  and  electromagnetic 
waves,  researchers  have  attempted  to  find  the  counterpart  left-handed  acoustic 
metamaterials,  namely,  materials  having  negative  mass  density  and  negative  elastic 
modulus  simultaneously.  Nevertheless,  in  contrast  to  the  natural  negative  electric 
permittivity,  no  natural  material  has  been  found  to  possess  either  negative  mass  density  or 
negative  modulus.  This  makes  the  realization  of  acoustic  counterpart  metamaterials  more 
difficult.  To  this  end,  some  researchers  have  resorted  to  man-made  microstructures  to 
achieve  an  apparent  effect  of  negative  mass  density  (Liu  et  al.,  2005;  Milton  and  Willis, 
2007;  Yao  et  al.,  2008;  Huang  et  al.,  2009;  Huang  and  Sun,  2009;  Lee  et  al.,  2009a). 
Most  of  these  researchers  employed  mechanical  resonators  as  “microstructures”  to  form 
the  acoustic  metamaterial,  which  was  first  introduced  by  Liu  et  al.  (2000).  It  was  found 
that  if  the  acoustic  metamaterial  is  treated  as  a  classical  elastic  solid,  then  its  effective 
mass  density  would  be  frequency-dependent  and  become  negative  in  certain  frequency 
range.  Recently,  acoustic  metamaterials,  possessing  negative  effective  modulus,  were 
proposed  and  studied  (Fang  et  al.,  2006;  Lee  et  al.,  2009b).  These  “metamaterials”  are, 
however,  structures  having  the  form  of  a  pipe  with  lateral  Helmholtz  resonators  or  side 
holes  designed  to  manipulate  acoustic  waves.  Attenuation  of  waves  was  also  confirmed 
experimentally  within  the  negative  range  of  effective  modulus. 
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It  is  also  of  interest  to  note  a  type  of  composite  materials  with  negative-stiffness 
inclusions  proposed  by  Lakes  et  al.  (2001)  in  which  the  working  mechanism  of  this  type 
of  materials  is  fundamentally  different  than  that  of  the  aforementioned  metamaterials.  In 
Lakes’  case,  negative  stiffness  can  occur,  for  instance,  when  the  deforming  object  has 
stored  energy.  This  property  is  usually  unstable  but  has  been  shown  theoretically  that 
negative  stiffness  can  be  stabilized  within  a  positive-stiffness  matrix  (Drugan,  2007). 
They  (Lakes  et  al.,  2001;  Drugan,  2007)  found  that  the  resulting  composites  exhibit 
extreme  mechanical  damping  and  large  anomalies  in  stiffness,  as  a  consequence  of  the 
high  local  strains  that  result  from  the  inclusions  which  deform  more  than  the  composite 
as  a  whole.  Jaglinski  et  al.  (2007)  further  showed  that  this  kind  of  composite  materials 
exhibit  a  viscoelastic  modulus  substantially  greater  than  that  of  diamond. 

In  the  present  paper  we  focus  on  the  acoustic  metamaterials  of  the  local-resonance 
kind.  First,  a  model  that  exhibits  negative  effective  Young’s  modulus  is  proposed  and  its 
dynamic  characteristics  are  investigated.  The  wave  propagation  behavior  of  the  one¬ 
dimensional  system  is  subsequently  studied.  Some  numerical  simulations  on  wave 
attenuation  and  filtering  in  the  acoustic  metamaterials  are  presented  at  the  end. 

4.2.  A  mechanical  model  with  extreme  and  negative  effective  Young’s  modulus 

4.2.1  A  unit  model  and  its  effective  Young’s  modulus 


Figure  4. 1  A  unit  of  the  mechanical  model  and  the  effective  continuum 


Consider  a  mechanical  model  shown  in  Figure  4.11  representing  a  unit  of  a  ID 
acoustic  metamaterial.  A  spring  with  elastic  constant  kj  and  length  L  in  its  equilibrium 
state  is  connected  to  four  rigid  and  massless  truss  members.  The  main  body  is  assembled 
symmetrically  as  shown  to  support  two  sets  of  microstructures  consisting  of  springs  with 
elastic  constant  k2  and  masses  m2.  Point  mass  m2  is  allowed  to  move  only  in  y  direction. 
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Consider  a  dynamic  test  of  stress-strain  relation  by  applying  a  symmetrical  pair  of 
harmonic  forces  F  in  the  x  direction  as  shown  in  Figure  4.1.  The  equations  of  motion  of 
the  system  are  readily  obtained  as 


m. 


at' 


+^^2(^2  -Vl) 

(4.1) 

■■k^{v,  -vj 

(4.2) 

In  the  equations  above,  the  external  force  assumes  harmonic  complex  form  F  -  Fe  ^ 

where  F  is  the  excitation  force  and  co  is  the  excitation  frequency.  Assuming  small 
displacements,  we  have 


(L-Ui)Ui  (L-Mj)'m' 


D 


■  + 


-I-...  =  -- 


D 


(4.3) 


The  force-displacement  relation  is  obtained  by  letting  the  displacement  fields  be 

_  ^  iOX  _  ^  iOX 

^1  “  ^1^  and  ^2  “  ^2^  and  substituting  the  harmonic  displacements  in  equations  (4.1) 
and  (4.2).  We  obtain 

,2 

’M, 

(4.4) 


F  =  2Ui  + 


1  k^co^ 

2  0}^  -  col 


If  the  mechanical  system  is  to  be  represented  by  a  ID  elastic  solid  as  shown  in 
Figure  1 ,  the  stress-strain  relation  of  the  effective  continuum  is 


^  ~  ^eff  ^ 


(4.5) 


^  -  F I A  (/i  is  the  cross-section  of  the  elastic  solid)  and  £  -u!  L  are  average  stress  and 

£ 

strain,  respectively,  of  the  continuum  element,  and  is  the  effective  Young’s  modulus 
of  the  ID  elastic  solid.  Since  the  models  in  the  present  study  are  ID,  in  the  following  the 
term  “modulus”  is  referred  to  “Young’s  modulus”  unless  specified  otherwise.  To  make 
the  two  models  dynamically  equivalent,  equations  (4.4)  and  (4.5)  must  be  the  same  and, 
thus,  the  effective  modulus  of  the  elastic  solid  would  take  the  form 


1  + 


d^i  rj^ 


^0  2/7  1  ^4 

where  EQ=k^LIA  is  the  static  Young’s  modulus,  ^  0)1  COq  non-dimensional  wave 

frequency  in  which  d  =  k^l  k^  ^  =  Df 

2.2  Discussion 

From  equation  (4.6),  we  note  that  the  effective  modulus  is  frequency-dependent  and 
is  unbounded  when  the  driving  frequency  ®  is  equal  to  the  local  resonance  frequency. 
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i.e.,  ^  .  The  dynamic  stress-strain  curves  illustrating  the  effective  moduli 

corresponding  to  different  frequencies  are  schematically  plotted  in  Figure  4.2.  In  the 
figure,  the  frequency-dependent  effective  modulus  is  seen  to  approach  infinity  when 

^  ! ^2  ^  and  approaches  zero  when  a*  ^vith 

of  =  (oj2/i2  +  Su) 

For  frequencies  in  between  coo  and  co  *,  the  effective  modulus  becomes  negative. 

Figure  3  shows  an  example  of  the  non-dimensionalized  effective  modulus  and  the 
displacement  ratio  V2/M1  with  respect  to  excitation  frequency.  In  the  example,  the  material 
and  geometrical  constants  are  selected  as:  =  ^2  ^ ,  /li  =  {LI D)^  =  1 .  For  harmonic 

vibrations,  the  displacement  ratio  V2/M1  is  obtained  from  equation  (4.2)  as 
V2  _ 

Ml  7^-1  (48) 

It  is  clearly  seen  that  for  frequencies  between  (Oo  and  co  *,  the  effective  modulus  and  the 
ratio  Vj/mi  would  become  negative. 


Figure  4.2  Stress-strain  curves  for  various  frequencies.  Go,  £0,  and  Eq  are  values  for  the 

corresponding  static  case 
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Figure  4.3  Dynamic  behavior  of  the  model  with  respect  to  excitation  frequency 

4.3  Significance  of  negative  effective  Young’s  modulus 

4.3.1  ID  lattice  model 


Consider  a  ID  model  consisting  of  a  large  number  of  resonator  units  as  shown  in 
Figure  4.4.  The  mass  of  the  background  material  in  this  ID  model  is  lumped  to  nodes  of 
the  lattice  system  and  denoted  by  m\.  The  equations  of  motion  for  the  j*  unit  are 


a2„o) 


m, 


dt^ 


—  kj  (m{-'  )  +  k^ 

—  Mp'  )  —  2^2  (v2^  ''  —  Vj  ■'  ''  )—— 

'  '2D 

-1-  2^2  {v2^  -  )-— 

2V  2  1 

(4.9) 

=  ^2  -  V2'''  ) 

(4.10) 

—u 

In  deriving  equations  (4.9-  4.10),  the  geometrical  relationship, 


is  used.  For  the  harmonic  wave  motion  given  by 


u 


(j+n)  _ 


=  u^e 


and 


V, 


(j+n)  _ 


we  derive  the  dispersion  relation  from  the  equations  of  motion 


with  the  result: 


S  2 


(l-coS(^)  \ri^  -1-2— (1- cos  (^)  =  0 


(4.11) 


where  d  =  m2lm^  ,  ^  is  non-dimensional  wave  number,  ^  -xlL  j 


is  non- 


dimensional  distance,  and  ^  is  non-dimensional  time.  The  dispersion  curves  are 

plotted  in  Figures  4.5  and  4.6,  respectively,  for  two  cases: 

Case  \:  e  =  2,  d  =  Q.5,  id  =  A. 
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Case  2:  ^  =  0.5,  ^  =  5,  //  =  1. 

In  either  case,  a  band  gap  associated  with  the  local  resonance  is  present.  Unlike  the  mass- 
in-mass  lattice  system  studied  by  Huang  et  al.  (2009),  the  local  resonance  frequency  of 
the  present  system  is  located  at  the  upper  bound  of  the  band  gap.  For  frequencies  in  the 

band  gap,  the  wave  number  must  assume  a  complex  form  ^  ^  in  order  to  satisfy 

the  dispersion  equation  (4.11).  This  implies  that  waves  with  a  frequency  in  the  band  gap 
would  attenuate  spatially.  It  is  noted  that,  in  the  present  case,  the  maximum  attenuation 
factor  p,  does  not  occur  at  the  local  resonance  frequency. 


It  is  of  interest  to  represent  the  ID  lattice  system  by  a  ID  elastic  solid  with  an 
effective  Young’s  modulus  Eeff-  Since  the  dispersion  curves  depict  the  dynamic  behavior 
of  the  ID  lattice  system,  the  representative  elastic  solid  should  match  the  dispersion 
behavior.  First,  take  the  effective  mass  of  the  system  to  be  the  average  mass  of  the  unit. 


namely,  classical  wave  equation  of  an  elastic  solid,  the 

dispersion  relation  is  given  by  .  Next,  substitute  this  dispersion  relation  in 

equation  (4.1 1)  to  obtain 


S  9  1 

-7^{cos-‘[l- 

u 


(77^-1)77^ 


2[{eid+eiiii2)T]^ -Old] 


]} 


-2 


(4.12) 


Alternatively,  we  can  adopt  the  effective  modulus  given  by  equation  (4.6)  which  is 
obtained  from  the  unit  (see  Figure  4.1)  by  a  dynamic  stiffness  equivalence  method. 
Effective  moduli  given  by  equations  (4.6)  and  (4.12)  are  plotted  in  Figures  4.5  and  4.6  for 
Casel  and  Case  2,  respectively.  The  two  methods  match  well  for  frequencies  up  to  a 
certain  range.  Since  the  dispersion-curve-matching  method  gives  the  exact  dispersion 
relation  of  the  lattice  model,  it  is  seen  that  the  dynamic  behavior  near  the  band  gap  can  be 
captured  accurately  by  the  corresponding  effective  continuum. 
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(Case  1) 


Figure  4.5  Results  for  case  1.  (a)  Dispersion  relation  (AC:  acoustic  mode;  OP:  optical 
mode);  (b)  attenuation  factor;  (c)  effective  modulus  obtained  by  (Ml)  dynamic  stiffness 
equivalence  method  of  a  unit,  and  by  (M2)  matching  of  dispersion  curves. 


(Case  2) 


Figure  4.6  Results  for  case  2.  (a)  Dispersion  relation  (AC:  acoustic  mode;  OP:  optical 
mode);  (b)  attenuation  factor;  (c)  effective  modulus  obtained  by  (Ml)  dynamic  stiffness 
equivalence  method  of  a  unit,  and  by  (M2)  matching  of  dispersion  curves. 
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4.3.2  Mechanism  of  maximum  wave  attenuation 


Of  particular  interest  is  the  mechanism  of  maximum  wave  attenuation  that  occurs 

at  the  frequency®  =  ®*.  Let  us  consider  a  unit  of  the  lattice  shown  in  Figure  4.7.  One 
side  of  the  unit  is  hinged  and  the  other  side  is  subjected  to  an  excitation  force  F.  The 
equations  of  motion  for  the  system  are: 

(mj  /  2)Mj  =  —  271  +  F  (4  13) 

m2  V2  =^2(^1 -^2)  (4.14) 

where  7;  is  the  horizontal  component  of  axial  force  7  in  the  massless  rigid  truss  member. 
From  Figure  7,  it  is  readily  seen  that 
F  =  Vi+27i 

7j  =  —72  =  —^2(^2  “^1) 

‘  D  2D 

in  which  T2  is  the  vertical  component  of  7. 


(4.16) 


Fj  — - Uy 

In  the  following  analysis,  the  geometrical  relationship,  2D  ^  which  is  valid 

for  small  displacement,  is  adopted.  For  steady  state  harmonic  motions,  F  =  Fe  ^ 

_  ^  iOX  _  ^  iOX 

“1  “  ^  and  ^2  -  ^2^  ^  and  from  equations  (4.13-  4.14)  we  obtain 


(4.17) 

The  time- 


—  m2(0^  +^2 
k^  L 

KJ 

1  A 

L 

.  Td  . 

A  1  4 

A  —  —m^m20)  — 


where 


1  1  212 

—  mjfcj +^2^1 +^'”2^2(2'/ f^)  +k^k2 


independent  reaction  force  F  is  obtained  by  substituting  equation  (4.17)  into  equations 
(4.15-  4.16)  and  using  equation  (4.7).  We  have 

R  _  m2k^{l  +  Sju/2)\ 


0)  +0)' 


■-] 


(4.18) 


It  is  noted  that  at  frequency  co  =  co*,  reaction  force  ^  =  0.  In  other 

words,  no  force  is  generated  at  the  far  side  of  the  unit.  This  explains  why  at  this 
frequency  this  lattice  yields  the  maximum  wave  attenuation. 
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Figure  4.7  A  single  unit  of  the  ID  lattice  model 
4.3.3  Possible  practical  designs 

One  question  arisen  is  how  to  practically  realize  acoustic  metamaterials  with 
negative  elastic  modulus.  In  Figure  8,  we  schematically  illustrate  one  of  the  possible 
designs  of  acoustic  metamaterials  with  frequency-dependent  elastic  modulus.  This  design 
take  directly  the  form  of  a  host  material  embedded  with  microstructures  of  the  kind 
shown  in  Figure  1.  The  spring,  ki  ,  now  represents  the  stiffness  of  the  base  material. 
Another  possibility  is  shown  in  Figure  9.  A  base  plate  is  employed  and  attached  with 
wavy  springs.  Concentrated  masses  are  placed  on  top  of  the  wavy  springs,  and,  as  a 
whole,  forming  the  lateral  local  resonators.  More  creative  designs  are  possible.  However, 
details  are  beyond  the  scope  of  the  present  study. 


<?>  4 


Figure  4.8  possible  design  of  metamaterial  with  frequency-dependent  effective 
modulus 


11 


Top  view 


Figure  4.9  2"*^  possible  design  of  metamaterial  with  frequency-dependent  effective 
modulus 

4.4  Dynamic  behavior  of  the  metamaterial 

4.4.1  Semi-infinite  ID  system 

Using  the  finite  element  method  we  perform  a  series  of  numerical  simulations  of 
sinusoidal  wave  propagation  in  a  ID  system  (see  Figure  4.4).  With  a  large  number  of 
lattice  units  to  prevent  the  reflection  of  waves,  the  system  can  be  considered  semi¬ 
infinite.  Numerical  results  are  obtained  using  the  commercial  finite  element  code 
ABAQUS  Explicit.  The  material  constants  of  Case  2  in  Section  4.3.1  are  chosen  so  that 

of  -  0.5345tyQ  members  are  assigned  with  a  high  stiffness  and  low  mass 

relative  to  the  springs  and  point  masses  so  that  they  can  be  regarded  as  rigid  truss 
members.  A  transient  harmonic  excitation  F  -  Fq  sin  ox,  t  >0,  is  applied  at  one  end  of 
the  lattice  system.  Various  input  frequencies  are  considered. 

4.4. 1.1  Extreme  effective  moduli 

It  is  of  interest  to  study  the  dynamic  stress-strain  relation  at  the  first  unit  when  the 
lattice  system  is  excited  by  an  external  force.  The  time  histories  of  the  stress-strain  curves 
are  plotted  in  Eigure  8.  Six  cases  with  different  excitation  frequencies  are  demonstrated. 

/  £ 

In  the  figure  the  horizontal  axis  denotes  non-dimensional  strain  °  of  the  first  unit, 

£ 

in  which  o  is  the  strain  when  force  Fq  is  applied  statically,  and 

where  is  the  strain  for  the  j*  unit,  and  is  the  displacement  of  mi  at  the  j*  unit.  In 
Eigure  4.10,  the  vertical  axis  is  the  non-dimensional  stress  produced  by  the  external  force. 
In  the  figure,  an  orange  dashed-line  denotes  the  entire  response  during  the  simulation 
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time  window  while  a  blue  solid  line  shows  the  last  few  cycles  of  the  response.  It  is 
evident  that  due  to  the  transient  effect,  the  orange  dashed  lines  appear  irregular  initially 
but  converges  to  the  steady-state  response  after  a  number  of  cycles  of  loading.  We 
particularly  note  the  dynamic  characteristics  of  the  first  lattice  unit  (j  =  1)  when  it 
exhibits  an  extreme  effective  modulus.  As  the  excitation  frequency  approaches  the  local 
resonance  frequency  coo,  the  dynamic  effective  modulus  of  the  unit  becomes  very  large, 
resulting  in  small  strains.  On  the  other  hand,  the  effective  modulus  becomes  small  as  the 
excitation  frequency  is  near  ^  * . 


(1)  T|  =  0.0316  (0  <  CO  <  «*) 


(5)  T|  =  0.9171  (oo  near  co^j) 
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(2)  T|  =  0.3162  (0  <  (0  <  CO*) 


Figure  10  Various  time  history  responses  of  the  stress-strain  relations  of  the  first  lattice 
unit  for  different  excitation  frequencies.  Orange  dashed-lines  represent  the  entire 
responses  of  the  simulation  time  window.  Blue  solid  lines  show  responses  for  the  last  few 
cycles  of  loading. 


4.4. 1.2  Spatial  wave  attenuation 


It  is  shown  in  Figures  5c  and  6c  that,  for  waves  with  frequencies  in  the  band  gap, 

the  dynamic  effective  modulus  becomes  negative.  To  further  understand  the  mechanics  of  negative 
modulus,  we  perform  numerical  simulations  for  frequencies  inside  and  outside  the  band  gap.  We  measure 
the  attenuation  of  wave  amplitude  by  considering  the  ratio  of  strains  in  the  j*  lattice  unit  and  (j+n)*  unit  as 


^(y) 

(4.20) 
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Figure  4.11  shows  the  results  for  the  first  six  lattice  units  for  various  frequencies.  The 
three  cases,  q  =  0.0316,  0.3162,  1.1068,  are  propagation  modes  while  q  =  0.5376,  0.7906, 
0.9171,  are  attenuation  modes.  It  is  clear  that  the  wave  amplitude  decays  spatially  when 
the  wave  frequency  is  in  the  band  gap,  especially  when  it  is  near  co*. 


Figure  4.11  Dimensionless  strain  envelopes  of  propagating  wave  fields  in  various 

frequencies,  coq  is  the  local  resonance  frequency;  co*  is  defined  in  equation  (7); 
represents  the  strain  of  the  n*  unit. 

4.4.2  A  metamaterial  sandwiched  between  two  elastic  media 


Medium  II; 
Metamaterial 


Figure  4.12  A  simple  illustration  showing  wave  propagation  through  a  metamaterial 
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sandwiched  between  two  elastic  media 


In  this  section  we  perform  numerical  simulations  of  transient  wave  propagation 
through  a  metamaterial  (Medium  II)  sandwiched  between  two  media  (Medium  I  and  III) 
as  depicted  in  Figure  4.12.  This  example  is  to  further  demonstrate  the  wave  blocking  and 
filtering  capability  of  the  present  metamaterial.  Using  a  representative  ID  mass-spring 
lattice  system  (see  Figure  4.4),  we  choose  twelve  units  for  the  metamaterial  in  the  x 
(longitudinal)  direction.  Each  unit  has  the  dimension  in  the  x  direction  equal  to  20  mm  so 
that  the  total  thickness  of  the  metamaterial  in  this  direction  is  0.24  m.  The  twelve  units 
are  divided  into  four  sections  with  different  units  of  properties.  In  other  words,  the  unit 
resonators  in  each  section  are  identical  but  differ  among  different  sections.  We  select  the 
material  parameters  and  corresponding  frequencies,  coq  and  co*,  of  each  section  as  listed 

in  Table  4.1  with  “  2.4e-5  kg  k  -  \.15e5  NIm  resulting  band  gap  of  the 
designed  metamaterial  is  covered  from  approximately  4100  Hz  to  7550  Hz. 

For  simplicity  and  without  loss  of  generality,  we  employ  a  simple  spring-mass 
lattice  system  to  represent  Medium  I  and  Medium  III.  In  the  simulation,  the  mass  and  the 

spring  of  this  representative  lattice  are  assumed  to  be  ^  kg 

k  -  l.leS  N / m  ^  respectively.  The  lattice  constant  is  ^  . 

Table  4.1  Material  constants  and  the  corresponding  frequencies  of 
microstructures  in  each  section  of  the  sandwiched  metamaterial 


Section 

1 

2 

3 

4 

0 

1.667 

2.000 

2.500 

3.000 

5 

0.514 

0.457 

0.429 

0.371 

F 

1.778 

1.778 

1.778 

1.778 

0)0  (Hz) 

7550 

6500 

5630 

4785 

0)*  (Hz) 

6255 

5480 

4790 

4145 

A  wave  packet  with  an  assigned  central  frequency  is  sent  from  the  far  side  to  the  left 
of  Medium  I  with  a  displacement-controlled  excitation.  The  central  frequency  of  the 
wave  packet  can  be  selected  so  that  its  frequency  spectrum  falls  inside  or  outside  the 
band  gap.  Two  cases  are  considered.  In  the  first  case,  the  wave  packet  with  the  central 
frequency  3185  Hz  is  generated.  The  frequency  band  starts  roughly  from  2500  Hz  to 
4100  Hz  so  that  it  does  not  contain  frequencies  in  the  band  gap  of  the  metamaterial. 
Snapshots  of  the  wave  packet  at  t  =  0.005  sec  and  t  =  0.050  sec  are  shown  in  Figure  4.13. 
D*  denotes  the  non-dimensionalized  amplitude  of  the  displacement  with  respect  to  the 
maximum  displacement  of  the  input  pulse.  It  is  obvious  that  the  wave  packet  propagates 
through  the  metamaterial  nearly  without  changing  in  shape  and  amplitude.  A  small 
amount  of  reflection  is  due  to  the  small  impedance  mismatch.  This  can  also  be  observed 
from  the  frequency  spectra  of  the  incident  and  transmitted  waves  shown  in  Figure  4.14. 

In  the  second  case,  the  wave  packet  with  the  central  frequency  5900  Hz  is  generated. 
The  frequency  band  is  roughly  from  4100  Hz  to  7500  Hz.  Thus,  it  covers  the  band  gap  of 
the  metamaterial.  The  snapshots  at  t  =  0.005  sec  and  t  =  0.050  sec  are  shown  in  Figure 
4.15.  It  is  evident  that  a  large  portion  of  the  wave  packet  is  reflected  by  the  metamaterial 
and  only  a  small  portion  of  the  wave  is  transmitted  into  the  Medium  III.  Moreover,  the 
frequency  spectra  of  the  incident  and  transmitted  waves  of  the  second  case  are  shown  in 
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Figure  4.16.  It  is  observed  that  a  large  amount  of  frequency  content  is  not  transmitted. 
Using  the  frequency  defined  in  equation  (7)  and  calculated  in  Table  4.1  for  each  section, 
we  denote  co*j  in  the  figure  where  “j”  is  the  section  number  of  the  sandwiched 

metamaterial.  It  is  seen  that  the  frequencies  with  maximum  wave  attenuation  match  well 
with  those  defined  by  equation  (4.7). 


(a)  Time  =  0.005s 
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(b)  Time  =  0.050s 
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Figure  4.13  Snapshots  of  the  wave  packet  with  a  central  frequency  of  3185  Hz  (outside 
the  band  gap  of  the  metamaterial)  (a)  before  and  (b)  after  it  propagates  through  the 
metamaterial.  The  metamaterial  is  located  at  x  =  40.00~40.24  m 


Figure  4.14  Frequency  spectra  of  the  wave  packet  with  a  central  frequency  of  3185  Hz  at 
certain  positions  in  the  fluid-like  media  for  (a)  incident  wave  and  (b)  transmitted  wave 
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(a)  Time  =  0.005s 
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(b)  Time  =  0.050s 
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Figure  4.15  Snapshots  of  the  wave  packet  with  a  central  frequency  of  5900  Hz  (a)  before 
and  (b)  after  it  propagates  through  the  metamaterial 


Figure  4.16  Frequency  spectra  of  the  wave  packet  with  a  central  frequency  of  5900  Hz  at 
certain  positions  in  the  fluid-like  media  for  (a)  incident  wave  and  (b)  transmitted  wave. 
Frequency  defined  in  Eq.  (7)  and  calculated  in  Table  1  for  each  section  is  denoted  as  a)*j , 
where  j  is  the  section  number  of  the  sandwiched  metamaterial. 


4.5  Conclusion 

We  have  proposed  a  “microstructural”  model  that  can  be  used  to  form  an  acoustic 
metamaterial  exhibiting  unusual  dynamic  mechanical  behavior.  It  was  shown  that  the 
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effective  Young’s  modulus  of  the  acoustic  metamaterial  is  frequency-dependent  and  can 
become  theoretically  unbounded  or  vanishingly  small  at  certain  frequencies.  Moreover,  in 
certain  frequency  range,  a  band  gap  is  created  in  which  the  effective  Young’s  modulus 
becomes  negative.  Using  a  ID  system  consisting  of  a  large  number  of  the  aforementioned 
microstructural  units,  we  demonstrated  through  numerical  simulations  these  novel 
dynamic  characteristics.  This  confirms  that  it  is  possible  to  create  metamaterials  by 
design  of  artificial  microsctructures  to  achieve  certain  unusual  mechanical  properties  that 
are  not  usually  present  in  natural  materials.  By  combining  the  present  model  with  the 
existing  acoustic  metamaterials  with  negative  effective  mass  density,  acoustic 
metamaterials  having  simultaneously  negative  effective  mass  density  and  negative 
effective  Young’s  modulus  can  be  realized. 

5.  Wave  Propagation  in  an  Acoustic  Metamaterial  Having 
Simultaneous  Negative  Mass  Density  and  Young’s  Modulus 

5.1  Introduction 

For  the  substances  having  simultaneously  negative  electric  permittivity  (s)  and 
magnetic  permeability  (p),  many  unusual  properties  like  reverse  Doppler  effect,  reverse 
Cherenkov  radiation,  and  negative  index  of  refraction,  may  arise  (Veslago,  1968).  The 
idea  has  been  that  when  both  s  and  p  are  negative,  the  metamaterial  sustain  waves  with 
group  velocity  opposite  to  the  phase  velocity,  and  the  waves  can  still  propagate  in  such 
media  without  attenuation. 

Similar  phenomena  may  apply  to  the  counterpart  acoustic  metamaterials  owing  to 
the  mathematical  analogy  between  acoustic  and  electromagnetic  waves.  Theoretical  or 
experimental  attempts  have  been  made  for  investigation  of  acoustic  metamaterials 
possessing,  for  instance,  negative  effective  mass  density  (Liu  et  al.,  2000;  Liu  et  al.; 
Milton  and  Willis,  2009;  Yao  et  al.,  2008;  Huang  et  al.,  2009;  Huang  and  Sun,  2009; 
Huang  and  Sun,  2011;  Lee  et  al.,  2009a,  Islam  and  Newas,  2012),  or  negative  effective 
elastic  modulus  Fang  et  al.,  2006;  Lee  et  al.,  2009b;  Huang  and  Sun,  2011),  or  both 
simultaneously  (Li  and  Chan,  2004;  Ding  et  al.  2007;  Cheng  et  al.,  2008;  Lee  et  al.,  2010; 
Fok  and  Zhang,  2011;  Hu  et  al.,  2011).  For  metamaterials  with  double  negativity,  the 
consequence  is,  as  mentioned,  metamaterial  behavior  of  negative  refraction,  and, 
equivalently,  negative  phase- velocity  and  backward-wave  phenomena  (Veslago,  1968; 
Lakhtakia,et  al.,  2002;  Lindell  et  al.,  2001;  Aydin  and  Guvan,  2005).  Nevertheless,  it  has 
been  shown  that  the  negative  refraction  can  be  achieved  not  only  by  metamaterials  with 
double  negative  properties,  but  also  by  Bragg  scattering  in  periodic  lattices  (Notomi, 
2000)  and  others  (Chen,  2008;  Ao  and  Chan,  2010).  The  similarity  between  these 
mechanisms  relies  on  negative  slope  of  dispersion  curves  in  the  first  quadrant  of 
dispersion  relation. 

According  to  the  definition,  the  negative  slope  of  the  dispersion  curves  has  been 
related  to  negative  group  velocity  in  anomalous  dispersion.  However,  as  we  agree  with 
some  references,  and  as  it  will  be  shown  in  the  present  study,  the  term  “negative  phase 
velocity”  is  better  describing  the  negative  slope  of  the  dispersion  curves  than  “negative 
group  velocity”,  at  least  for  the  acoustic  metamaterials  with  double  negativity.  It  also 
distinguishes  the  aforementioned  metamaterials  from  the  negative  group-velocity  media 
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resulting  in  fast-light  phenomena  (Wang  et  ah,  2000;  Boyd  and  Gauthier,  2002; 
Feigenbaum,  2009;  Boyd  and  Gauthier,  2009).  It  is  noted  that  the  subject  of  the  negative 
phase  velocity  has,  surprisingly,  been  discussed  since,  at  least,  1904  by  Lamb  (Lamb, 
1904),  using  mechanical  models,  and  by  Schuster  (Schuster,  1904)  in  electromagnetism. 
At  that  time,  the  researchers  were  pessimistic  about  the  physical  applications. 

The  study  of  EM  waves  in  anomalous  dispersion  has  been  extensively  carried  out  by 
researchers  in  recent  years,  but  that  of  acoustic  and  elastic  waves  receives  much  less 
attention.  In  the  present  study  we  propose  a  one-dimensional  (ID)  elastic  mass-spring 
system  consisting  of  material  units  exhibiting  simultaneously  negative  effective  mass 
density  and  effective  elastic  modulus  in  a  spectral  band.  The  local  resonance  acts  as  the 
central  mechanisms  of  the  unusual  behavior.  Special  attention  is  focused  on  the 
frequency  band  exhibiting  negative  slope  of  the  dispersion  curves.  Transient  wave 
propagation  is  investigated  computationally  to  better  understand  the  dynamic 
characteristics  regarding  negative  phase  velocity. 

5.2  Theoretical  models  of  acoustic  metamaterials 

5.2.1  Model  of  negative  mass  density  (NMD) 

Let  us  begin  with  the  mass-in-mass  lattice  model  that  has  been  discussed  by 
Huang  and  Sun  (2009).  As  illustrated  in  Figure  la,  the  mass-in-mass  unit  takes  the  form 
of  a  rigid  ring  with  mass  nij  which  is  connected  to  the  neighboring  units  by  springs  with 
spring  constant  kj.  The  ring  contains  an  internal  mass  m2  connected  to  the  ring  by  an 
internal  spring  with  spring  constant  k2.  Both  m;  and  m2  are  only  allowed  to  move  in  the  x 
(horizontal)  direction. 


(a) 


^NMD 

St 


y 


(b) 


Figure  5.1  (a)  One-dimensional  lattice  model  and  the  continuum  representative  having 
negative  effective  mass  density  (NMD);  and  (b)  One-dimensional  lattice  model  and  the 
continuum  representative  having  negative  effective  Young’s  modulus  (NYM). 


The  equation  of  motion  for  the  j*  unit  cell  are  given  by 
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dt^  IV  1  1 


)+^2(mi^'  -uY^)=0 


h - ^ 

dt 


where  ^  (y  =  1  and  2)  denotes  the  displacement  of  mass  ‘y’  in  the  j*  cell.  The 
harmonic  wave  solution  for  the  (j+n)*  unit  is  expressed  by 

,.(i+n)  _  ,7  Ji^'  +n4-rfT‘‘ ) 

-u^e  ^^2) 

In  Eq.  (3),  ^  ~  is  the  non-dimensional  wavenumber,  rj  -  (ol (0^^°  the  non-dimensional 


_  i^C  /  TYl  _ 

wave  frequency,  *"0  ^  2  the  local  resonance  frequency,  ^  =xil  jhe  non- 

Y*  —  f 

dimensional  spatial  parameter  with  the  lattice  spacing  L,  and  ~  o  the  non- 
dimensional  time.  The  dispersion  relation  is  obtained  by  substituting  the  harmonic 
displacements  in  Eq.  (5.3)  and  solving  the  resulting  eigen-value  problem  of  the 
coefficients.  We  obtain 

=  ri^  -  (l  +  ^2i)  +  ^^(l-cos,^)  +^^(l-cos,^)  =  0 


the  local  resonance  frequency,  ^  -xlL  the 


where  6^^  =mj^lm^  and  <^21  =k2lk^  are  the  non-dimensional  mass  ratio  and  stiffness 
ratio,  respectively. 

An  effective  elastic  solid  was  developed  by  Huang  and  Sun  (2011b)  to  represent  the 
original  mass-in-mass  lattice.  It  was  found  that  the  effective  mass  density  of  the 
equivalent  elastic  solid  is  frequency-dependent  in  the  form 

(  0  77^  ^ 

NMD  _  NMD  I  .  '^21  '/ 

Peff  -Pst  1_„2 

V  ^  +  ^21  1  V  J 

where  =  (^n^  +  m2) /  AL  is  the  static  mass  density  in  which  A  and  L  are  the  cross- 

sectional  area  and  the  lattice  spacing,  respectively.  This  model  is  denoted  as  “NMD” 
since  it  exhibits  negative  effective  mass  density  in  the  range 


The  time-independent  wave  equation  for  plane  harmonic  waves  is  simply  given  by 

2  ^ 

2-  O  U 

Peff  U  —  ^ 

(5.7) 

p  = 

where  the  effective  mass  density  is  given  in  Eq.  (5.5),  and  the  effective 

E  =  =  k  L/  A 

Young’s  modulus,  ,  is  determined  from  the  static  stress-strain 

relation. 
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The  harmonic  wave  solution  for  the  equivalent  elastic  solid  is  given  by  u-ue 
Hence  the  dispersion  relation  is  obtained  as 


Gnmd ex  -NmdI-.nmdV 
elast  Peff  H  ) 


if  - 


^^11}  =0 


(5.8) 


5.2.2  Model  of  negative  Young’s  modulus  (NYM) 

Next  we  consider  the  mechanical  lattice  model  having  lateral  local  resonators 
shown  in  Figure  5.1b.  This  model  was  investigated  and  discussed  by  Huang  and  Sun 
(2011a).  Briefly,  point  masses  mi  are  connected  by  springs  with  a  spring  constant  ki. 
Four  rigid  and  massless  truss  members  are  assembled  symmetrically  as  shown  to  support 
two  sets  of  microstructures  consisting  of  springs  with  spring  constant  ks  and  masses  m^. 
Point  mass  ms  is  allowed  to  move  only  in  the  y  (vertical)  direction.  L  and  D  are 
geometrical  parameters  shown  in  the  figure. 

For  this  lattice  model,  the  equations  of  motion  for  the  j*  unit  are  given  by 


(5.9) 


m,^^  +  kM'^ -vf)  =  0 

dt  (5.10) 

From  geometrical  relations  based  on  the  assumption  of  small  displacements,  we  have 
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In  a  similar  manner,  the  dispersion  relation  is  obtained  as 


(5.11) 
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where  ^3^  =  Im^  ,  =k^/k^,  Tj  =  col  0)^  , 

The  relation  of  the  non-dimensional  wave  frequencies  reads 

7  =W(^31^2l)/(^2l‘^3l) 


= Vm 


(5.12) 


,  and  jii  =  {LIDf  . 


(5.13) 


If  this  lattice  system  is  represented  by  an  equivalent  elastic  solid,  the  effective 
Young’s  modulus  of  the  solid  is  found  frequency-dependent  in  the  form 
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^NYM 

where  '  is  the  static  Young’s  modulus.  This  model,  denoted  as  “NYM”, 

exhibits  negative  effective  Young’s  modulus  in  the  range 


21 


CO 


2  + 


<71^ 


<1 


CO' 


(5.15) 


The  wave  equation  takes  the  form  of  Eq.  (5.7),  in  which  the  effective  Young’s 


modulus, 


^eff  -  ^ 


NYM 

eff 


,  is  given  by  Eq.  (5.14)  and  the  effective  mass  density  is  the  static 


P.ff=pr=mJAL 


formulation. 


.  (Note  that  the  masses  mj  in  this  model  are  assumed 


to  move  only  in  the  y  direction  so  that  they  contribute  zero  inertia  forces  in  the  x 
direction).  The  dispersion  relation  is  given  by 

GNYM  -NYM  (  ..NYM  Y  jpNYM  Cl  i  j2  a 

elastiV^G)^Ps,  H  )V  -E,ff  ^  IL  ^0 


5.2.3  Model  of  double  negativity  (DN) 


Figure  5.2  One-dimensional  lattice  model  exhibiting  double  negativity  (DN) 


Consider  a  ID  infinite  lattice  shown  in  Figure  5.2.  This  model  is  a  combination 
of  the  two  models  in  Figures  5.1  such  that  mass  m;  in  the  NYM  model  is  replaced  by  the 
mass  unit  in  NMD.  This  model  is  denoted  as  “DN”  as  it  will  be  shown  to  possess  double 


negativity.  Let  the  displacements  of  masses  mi,  m2,  and  m3  at  the  j*  unit  be 
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^2  Ku  ^  3  \  2  ^  respectively.  The  equations  of  motion  of  the  system  at  the  j  cell  are 


readily  obtained  as 
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The  steady-state  solution  of  the  equations  of  motion  (5.17-5.19)  is  determined  by 
letting  the  displaeements  be 
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and  ‘^3  -  (5.20) 

where,  again,  ^  ~  is  the  non-dimensional  wavenumber,  ^  =  x! L  ^^d  ~  ® 

are  the  non-dimensional  parameters  in  space  and  time,  respectively.  Substituting  these 
harmonic  displacements  in  Eqs.  (5.17)  to  (5.19)  and  solving  the  resulting  eigen-value 
problem  of  the  coefficients,  we  obtain  the  dispersion  relation: 
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(5.21) 


On  the  other  hand,  the  equation  of  motion  of  the  equivalent  one-dimensional  elastic 
solid  is  given  by  Eq.  (5.7).  In  this  formulation  the  effective  mass  density  and  the 


effective  Young’s  modulus  are  given  by  Eq.  (5.5), 


^  ^NMD 

Heff  Hejf  Hejf 


and  Eq.  (5.14), 


ejf  eff  eff  ^  rcspectively.  Consequently,  the  dispersion  relation  of  the 
representative  elastic  solid  is  obtained  as 

GDN  ex  „DN  i  TpDN  t:2  !  j2  A 

elaAV4)^  }V  “  ^  =0  (5  22) 


5.3.  Effective  dynamic  parameters  and  dispersion  relation 

For  numerical  illustrations,  we  present  two  cases  with  the  following  assumed 
material  constants: 


(Case  1)  ^2i='«2/'«i=1-5,  ^ 
fi  =  {LIDf  =2.1% 


31 


m 


lm,=\  (?2i  =^2/^1  =0-1  <^31  =^3/^1=! 


and 


(Case  2) 


^21  =  2.25  ^31  =  0.4  ^21  =  0.3  =  0.25  //  =  25 


As  has  been  mentioned  in  Section  5.2.3,  the  model  in  Figure  5.2  is  denoted  as  “DN” 
because  for  some  range  of  wave  frequencies  it  could  exhibit  double  negativity.  This 
range  can  be  obtained  from  the  negative  modulus  and  mass  density  properties  and  their 
corresponding  spectral  regions.  The  effective  mass  density  in  Eq.  (5.5)  and  the  effective 
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Young’s  modulus  in  Eq.  (5.14)  with  respect  to  the  non-dimensional  wave  frequency  for 
Case  1  and  Case  2  are  illustrated  in  Figures  3  and  4,  respectively. 

In  Figure  5.3  for  Case  1,  the  material  constants  are  selected  so  that  the  frequency 
range  of  the  negative  quantity  of  the  NMD  model  lies  apart  from  that  of  the  NYM  model. 
In  other  words,  for  the  frequency  range  in-between  the  two  regions  with  negative  quantity 
(e.g.,  r\  approximately  ranges  in  1.6  and  2.5  in  Figure  5.3),  both  the  effective  mass 
density  and  the  effective  Young’s  modulus  are  positive.  On  the  other  hand,  in  Figure  5.4 
for  Case  2,  the  material  constants  are  chosen  so  that  the  frequency  ranges  of  both 
negative  quantities  of  the  NMD  and  the  NYM  models  overlap. 


Case  1 


Figure  5.3  Effective  dynamic  non-dimensional  parameters  for  Case  1  using  Eqs.  (5.5)  and 
(5.14).  The  bounds  of  negative  quantities  are  described  in  Eqs.  (5.6)  and  (5.15). 


Case  2 


Figure  5.4  Effective  dynamic  non-dimensional  parameters  for  Case  2.  Descriptions  are 
referred  to  the  caption  of  Figure  5.3. 
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Next,  the  dispersion  curves  for  both  cases  obtained  using  the  exact  formulation  (Eqs. 
(5.4),  (5.12),  and  (5.21)  for  NMD,  NYM,  and  DN  models,  respectively)  of  the  lattice 
model  are  shown  in  Figures  5.5  and  5.6,  respectively,  with  solid  lines.  In  each  figure, 
three  columns  labeling  DN,  NMD,  and  NYM  represent  the  dispersion  curves  of  the 
double  negativity  model,  the  negative  mass  density  model,  and  the  negative  Young’s 
modulus  model,  respectively.  For  each  diagram,  band  gaps,  where  waves  attenuated  and 
stopped,  exist  due  to  the  effect  of  the  local  resonance. 


DN  NMD  NYM 


Figure  5.5  Dispersion  curves  for  Case  1.  DN:  double  negativity  model;  NMD:  negative 
mass  density  model;  NYM:  negative  Young’s  modulus  model;  Solid  lines:  mass-spring 
lattice  model  from  Eqs.  (4),  (12),  and  (21)  for  NMD,  NYM,  and  DN  models,  respectively; 
Circles:  equivalent  elastic  solid  from  Eqs.  (8),  (16),  and  (22)  for  NMD,  NYM,  and  DN 
models,  respectively. 


The  accuracy  of  the  elastic  solid  representation  can  be  verified  by  comparing  the 
dispersion  curves  obtained  with  the  lattice  and  the  elastic  solid  models.  In  Figures  5.5 
and  5.6,  the  dispersion  curves  of  the  equivalent  elastic  solid  are  obtained  using  Eqs.  (5.8), 
(5.16),  and  (5.22)  for  the  NMD,  NYM,  and  DN  models,  respectively.  Good 
approximation  is  observed  for  the  lower  branch  for  all  models.  For  the  higher  branch 
some  deviation  is  seen  at  the  high  frequency  range.  The  reason  is  that  the  equivalent 
elastic  solid  employs  an  approximate  formulation  from  the  original  lattice  system  based 
on  a  long  wavelength  approximation.  It  is  noted  that  the  shaded  areas  in  the  figures 
indicate  the  negative  properties  derived  from  Eqs.  (5.6)  and  (5.15)  for  the  NMD  and 
NYM  model,  respectively.  For  both  cases  and  for  the  DN  model,  the  band  gaps  closely 
match  the  shaded  areas  of  the  separate  models. 


25 


DN  NMD  NYM 


Figure  5.6  Dispersion  curves  for  Case  2;  Solid  lines:  mass-spring  lattice  model;  Circles: 
elastic  solid  model;  detailed  description  referred  to  Figure  5.5. 


5.4  Phase  velocity  and  group  velocity 
5.4.1  Phase  velocity 


In  the  present  study  we  adopt  the  definition  of  the  phase  velocity  of  a  harmonic 
wave  as  the  following  form 

^  =- 

^  (5.23) 

All  points  of  the  same  phase  propagate  in  space  with  this  velocity.  In  the  dispersion  co-q 
diagram,  the  phase  velocity  is  given  by  the  slope  of  the  straight  line  connecting  the  origin 
and  a  certain  point  on  the  dispersion  curve.  Since  negative  wave  frequency  is  not  in  our 
consideration,  it  is  obvious  that  the  wave  vector,  q,  points  in  the  direction  of  phase 
velocity. 

For  an  effective  elastic  solid  of  an  acoustic  metamaterial,  the  dispersion  relation  can 
be  obtained  by,  for  instance,  employing  Eq.  (5.7)  together  with  in  Eq.  (5.5)  and 

in  Eq.  (5.14),  and  substituting  the  harmonic  plane-wave  solution  into  the  Eq.  (5.7).  The 
dispersion  relation  reads 

G{CO,q)=p^^CO^-E^ffq^=0 

From  the  dispersion  relation,  one  can  immediately  conclude  that  in  order  for  the  elastic 
medium  to  have  wave  propagation  without  attenuation  (i.e.,  a  real  solution  for  q),  the 
medium  must  have  either  both  positive  and  or  both  negative  and  E^g  . 
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The  phase  velocity  is  thus  given  by 


=  — =  + 


^ejf 


Pejf 


(5.25) 

This  is  still  not  possible  to  determine  the  direction  of  a  wave  in  the  acoustic  metamaterial 
from  the  sign  of  phase  velocity,  or  the  wave  vector,  based  on  the  given  dispersion 
relation,  since  the  real  wave  vector  q  could  be  either  positive  or  negative.  We  shall 
proceed  to  the  next  subsections  for  further  investigation. 


5.4.2  Group  velocity  of  DN  acoustic  metamaterial 


Of  particular  interest  is  to  investigate  the  sign  of  group  velocity  of  a  wave  group  or 
packet  in  an  acoustic  metamaterial.  This  wave  group  or  packet  corresponds  to  a 
superposition  of  harmonic  waves  with  very  similar  wave  frequencies  co  and  wavenumbers 
q.  In  the  dispersion  co-q  diagram,  the  group  velocity  of  a  wave  group  or  packet  is  given 
by  the  slope  of  the  dispersion  curve.  By  using  the  implicit  differentiation,  the  final 
expression  of  the  group  velocity  reads 
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From  the  expressions  given  by  Eq.  (5.5)  and  Eq.  (5.14)  it  can  be  shown  that 
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(5.28) 


where  ^  “  ^21  /(I  +  ^11)  ^  xhe  term 
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of  in  the  denominator  in  Eq.  (5.26)  is 
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always  positive  if  p^^  and  are  of  the  same  sign.  Substituting  Eqs.  (5.5)  and  (5.26)  in 
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V  / 

Since  0  <  /I  <  1  is  always  true,  this  inequality  holds  for  all  frequencies.  It  is,  therefore, 
proved  that  the  denominator  of  Eq.  (26)  is  positive  for  identical  signs  of  p^g  and  . 

Obviously,  from  Eq.  (5.26),  group  velocity,  ,  is  proportional  to  .  The 

sign  of  group  velocity  is,  hence,  affected  by  the  sign  of  the  effective  Young’s  modulus 
E^g  and  the  wave  number  q .  From  the  last  section,  in  order  for  an  elastic  medium  to 
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have  wave  propagation  without  attenuation,  the  medium  must  possess  either  both  positive 
p^g  and  E^g  or  both  negative  p^g  and  E^g  .  If  p^g  and  E^g  of  the  medium  are  both 

positive,  group  velocity  is  obviously  positive  for  a  positive  wavenumber.  On  the 
other  hand,  if p^g  and  E^g  are  simultaneously  negative,  waves,  therefore,  propagate  in 
the  medium  with  opposite  signs  of  phase  velocity  and  group  velocity  . 

5.4.3  Negative  phase  velocity 

For  waves  in  a  lossless  elastic  medium  the  energy  velocity  is  identical  to  the  group 
velocity  (Brillouin,  1953).  In  view  of  positive  energy  flow  for  a  lossless  elastic  medium, 
the  group  velocity  is  always  positive.  Once  again,  we  exclude  the  negative  group- 
velocity  media  resulting  in  fast-light  phenomena.  That  negative  group-velocity 
phenomenon  is  beyond  the  scope  of  the  present  study. 

For  the  acoustic  metamaterial  with  frequency-dependent  effective  mass  density  p^g 

and  Young’s  modulus  E^g ,  if  p^g  and  E^g  are  both  positive,  a  positive  group  velocity 
results  in  a  positive  wave  number.  On  the  other  hand,  if  Pg^  and  Eg^  are  negative 

simultaneously,  in  order  to  have  a  positive  group  velocity,  wave  number  must  be 
negative.  This  means  that  in  the  double  negativity  frequency  region,  the  sign  of  phase 
velocity  in  Eq.  (5.25)  must  be  negative. 

5.5  Dynamic  characteristics  of  wave  propagation  by  simulation 
5.5.1  Harmonic  waves 

To  investigate  how  a  group  of  waves  of  essentially  the  same  wavelength 
propagates  with  a  negative  phase  velocity,  we  perform  numerical  simulations  of  transient 
wave  propagations.  One  of  the  reasons  the  transient  wave  propagation  rather  than  steady- 
state  wave  propagation  is  studied  is  because  the  direction  of  disturbances,  hence  the 
direction  of  group  velocity,  can  be  ensured. 

Dynamic  simulations  are  implemented  by  use  of  the  ABAQUS  Explicit  solver.  In 
each  simulation,  a  lattice  model  shown  in  Figure  5.2  is  constructed  with  a  sinusoidal 

displacement,  sin(tyt)/f(0  ^  l^j^g  lattice.  The  material 

properties  of  Case  2  listed  in  Section  3  together  with  - 115.47 rad/ 5  used.  If 
the  lattice  system  is  long  enough,  the  reflected  wave  does  not  appear  in  the  window  of 
numerical  results. 
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Figure  5.7  Snapshots  of  sinusoidal  wave  propagation  in  finite  DN  model  with  excitation 
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Figure  5.8  Snapshots  of  sinusoidal  wave  propagation  in  finite  DN  model  with  excitation 
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Figures  5.7  and  5.8  show  the  snapshots  of  the  displacement  fields  in  space  with  two 

„  ^  j  ^NMD  ^ 

different  wave  frequencies,  ^  ~  ^  =  0.346  and  1.559,  respectively.  denotes 

the  normalized  amplitude  of  the  displacements  with  respect  to 
typical  propagation  mode  is  demonstrated.  Three  consecutive 


r  =tsAr 


with  the  non-dimensional  time  increment 


AT  =0.462 


In  each  figure,  a 
sections  of  time, 

show  the  behavior 
AT* 

of  wave  propagation.  For  now,  let  us  assume  that  the  time  increment  is  sufficiently 
small  in  capturing  the  dynamic  behavior  of  the  model.  We  shall  come  to  this  point 


m 


NMD 


t 


immediately  in  the  following  paragraph.  Note  that,  again,  ~  o  is  the  non- 
dimensional  parameter  in  time  domain  and  t  is  in  seconds,  “ts”  denotes  the  time  steps 
hereafter.  Section  1  (ts  =  1~4),  for  instance,  is  the  initial  state  and  the  section  3  (ts  = 
101-104)  is  considered  close  to  steady  state. 

AT* 

Selection  of  the  time  increment,  ,  for  visualization  of  simulation  results  is 
nontrivial.  An  improper  selection  of  the  increments  may,  on  the  one  hand,  produce 
incorrect  simulation  results,  or,  on  the  other  hand,  simply  cause  expensive  computation. 
Here  we  provide  a  simple  guide  for  the  balance  of  accuracy  and  efficiency.  The 
definition  of  the  time  increment  reads 


AT'*  ‘  ..Nl 

AT  =  COq  —  =  tyQ 

ts 


•At 


(5.30) 


where  At  is  the  time  step  in  seconds.  Substituting  the  non-dimensional  wave  frequency 

^  ^  j  ^NMD 

‘  0  in  Eq.  (5.30)  and  use  the  relationship  of  angular  frequency  and  period  of 


time  we  obtain 


(5.31) 


For  a  harmonic  wave.  At  <  /8  should  be  sufficient  in  capturing  the  dynamic  behavior 

of  the  models  in  the  present  study.  For  Case  2,  for  instance,  this  inequality  At  <  /  8 

and  the  wave  frequency  7-1-559  ^  together  with  Eq.  (5.31),  result  in  AT*  <0.504. 


Clearly  the  time  increments  chosen  for  both  cases  are  satisfactory. 

In  Figure  5.7,  given  the  first  section  from  ts  =  1  to  4,  the  propagation  direction  of  the 
front  of  the  dynamic  disturbance  is  seen  pointing  from  left  to  right,  indicating  the 
direction  of  the  energy  flow.  However,  due  to  the  nature  of  transient  wave  propagation, 
the  forerunner  of  the  dynamic  disturbance  may  contain  frequencies  other  than  the 
excitation  frequency  co.  Further  time  steps,  in  the  middle  and  bottom  windows  of  Figure 
5.7,  need  to  be  carried  out  in  order  to  make  sure  the  energy  direction  of  the  sinusoidal 
wave  with  exactly  the  excitation  frequency  (o.  In  the  second  section,  ts  =  51  to  54,  it  is 
seen  that  the  transient  part  of  the  wave  already  propagates  to  the  right  side  of  the  spatial 
window,  leaving  the  transition  to  steady  state  at  the  left  part  of  the  window.  From  the 
third  section,  ts  =  101  to  104,  a  more  complete  configuration  approaching  the  steady- state 
wave  propagation  can  be  observed.  At  this  stage,  the  propagation  direction  of  the  wave 
phase  is  observed  pointing  from  left  to  right. 
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Figure  5.8  shows  the  profile  of  wave  propagation  with  a  frequency  in  the  region  of 
double  negativity.  Noticeable  difference  to  Figure  5.7  is  the  negative  propagation 
direction  of  the  wave  phase  with  respect  to  positive  energy  flow.  In  the  second  and  third 
windows,  the  phase  direction  is  clearly  seen  pointing  from  right  to  left  during  the 
sequence  of  time  evolution. 

5.5.2  Wave  packets 


Group  velocity  is  the  velocity  of  the  envelope  of  a  wave  group  or  packet.  For  a 
wave  packet  corresponding  to  a  superposition  of  harmonic  waves  with  a  finite  frequency 
band,  the  propagation  behavior  in  a  dispersive  media  like  DN  acoustic  metamaterial 
should  be  of  particular  interest.  One  great  advantage  of  considering  the  wave  packet  is 
that  the  frequency  content  of  the  packet  can  be  controlled  and  the  energy  propagation 
direction  can  be  easily  identified  as  the  propagation  direction  of  the  wave  packet. 

To  this  end,  wave  packets  propagating  in  the  metamaterial  with  double  negativity 
are  simulated  in  a  similar  fashion.  Two  wave  packets  in  time  as  well  as  in  frequency 
domain  used  for  the  simulation  are  illustrated  in  Figure  9.  In  Figure  9a,  the  wave  packet 

possesses  a  frequency  band  approximately  from  V  —  —  0-22  7  -  0-47 

double  positivity  region)  with  the  central  frequency  at  ^  ~  0.346  ^  pjgure  9b,  the 

wave  packet  possesses  a  frequency  band  approximately  from  ^  ^^^o.nmd  1-10 

77-2.10  |.j^g  negativity  region)  with  the  central  frequency  at  ^  ~  1.559  p^ch 

packet  is  sent  from  the  left  end  of  the  model.  Same  material  properties  are  chosen  so  that 
the  dispersion  relation  of  the  model  is  shown  in  Figure  5.6. 


(a)  T*  =  ts  •  AT* 
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Figure  5.9  Two  typical  wave  packets  in  time  and  frequency  domain  with  (a)  frequency 
band  approximately  from  -0.22  ;7_o.47  jj^g  ^jo^p^g  positivity 

region)  with  the  central  frequency  at  ^  “  0.346  frequency  band  approximately 

from  V  -  ^0  -1.10  77-2.10  ^j^g  negativity  region)  with  the  central 

frequency  at  ^  ~  1.559 
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Figure  5.10  Snapshots  of  wave-packet  propagation  in  finite  DN  model;  Wave  packet 
consists  of  frequency  band  approximately  from  V  -  -  0*22  -  0.47 

central  frequency  at  ^  =  0-346  ^115  AT  rad  /  s 

region;  ts  =  time  step; 


(in  the  double-positivity 


1  =  (On  t  =  tS  -  A1 


AT  =0.462, 
with  ). 


Figure  5.10  shows  a  few  snapshots  of  wave  propagation  of  the  wave  packet  shown 
in  Figure  5.9a,  of  which  the  central  frequency  is  located  at  the  positive  effective  mass 
density  and  positive  effective  elastic  modulus  region  as  shown  in  Figure  5.4.  The  wave 
packet  propagates  without  much  distortion  since  most  of  the  spectral  content  travel  with 
approximately  the  same  speed,  which  is  also  evident  in  view  of  the  dispersion  curves 
shown  in  Figure  5.6.  The  energy  flow  direction  is  the  moving  direction  of  the  wave 
packet,  which  is  clearly  from  left  to  right.  Note  that  a  wave  packet  contains  not  just  a 
single  frequency,  but  a  finite  frequency  band.  It  is  therefore  difficult  to  identify  the  phase 
velocity  of  from  the  propagation  of  the  wave  packet. 

Next,  consider  the  case  when  the  frequency  components  of  the  wave  packet  are 
selected  in  the  frequency  region  of  double  negativity.  A  wave  packet  of  Figure  9b  is 
generated  for  simulation.  The  frequency  band  of  this  wave  packet  is  inside  the  double 
negativity  region  shown  in  Figure  5.4.  The  wave  packet  is  sent  from  the  left  end  and  is 
clearly  seen  to  propagate  in  the  positive  direction  (left  to  right)  from  the  top  and  bottom 
windows  of  Figure  5.11.  Negative  phase  direction  of  the  packet  is  expected  and  can  be 
verified  from  the  middle  window  of  Figure  5.11. 
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Figure  5.11  Snapshots  of  wave-packet  propagation  in  the  DN  model  Wave  packet  consists 
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5.5.3  Discussion 

The  region  of  the  dispersion  curve  with  negative  slope  means  a  spectral  range  where 
the  phase  and  group  velocities  are  opposite  in  direction.  Some  researchers  in  the 
literature  have  termed  the  backward-wave  phenomenon  along  with  “negative  group 
velocity”  simply  because,  perhaps,  in  the  co-q  diagram  the  slope  of  the  dispersion  curve  is 
negative.  For  example,  let  us  refer  to  the  dispersion  diagram  in  Figure  12  for  illustration. 
The  dispersion  curves  are  reproduced  from  Figure  6,  Case  2  in  Section  5,  for  the  DN 
model.  Usually  the  first  quadrant  of  the  diagram  is  considered  where  oi  and  q  are  both 
positive.  For  the  first  lowest  propagation  mode,  both  phase  velocity  and  group  velocity, 
according  to  the  definition  in  Section  4,  are  positive.  As  for  the  second  lowest 
propagation  mode,  in  the  first  quadrant,  the  phase  velocity  is  positive  while  the  group 
velocity  is  negative. 

The  misleading  and  confusing  conclusion  of  the  “negative  group  velocity”  may 
come  from  the  fact  that  energy  velocity  is  not  directly  provided  by  the  dispersion  curves 
in  the  first  quadrant.  In  fact,  the  second  branch  of  the  dispersion  curve  should  be  plotted 
in  the  second  quadrant  where  o)  is  positive  and  q  is  negative.  By  so  doing,  the  group 
velocity  is  positive  while  the  phase  velocity  is,  hence,  negative. 
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5.6  Conclusion 


We  have  proposed  a  one-dimensional  mechanical  model  representing  an  acoustic 
metamaterial.  The  system  can  be  modeled  with  an  equivalent  elastic  solid  which  exhibits 
simultaneously  negative  effective  mass  density  and  effective  Young’s  modulus  in  a 
certain  frequency  range.  In  this  double  negativity  frequency  range,  harmonic  waves 
propagate  without  attenuation.  However,  the  phase  velocity  becomes  opposite  to  the 
energy  flow  direction.  In  view  of  positive  energy  direction,  it  is  the  backward  wave  that 
propagates  with  negative  phase  velocity.  Transient  wave  simulation  results  not  only 
agree  well  with  the  analytical  prediction,  but  also  clearly  illustrate  the  unusual  wave 
phenomenon. 

6.  DESIGN  AND  ANALYSIS  OF  ANISOTROPIC  ELASTIC  METAMATERIALS 

6.1  Introduction 

An  anisotropic  elastic  metamaterial  made  of  heavy  cylinder  cores  coated  with 
elliptical  soft  layer  and  embedded  in  a  matrix  as  shown  in  Fig.  la  is  considered.  The 
representative  volume  element  (RVE)  is  identified  and  shown  in  Fig.  1(b).  The 
microstructure  is  distributed  in  a  rectangular  lattice  array.  In  the  figure,  the  hard 
inclusion  core  is  labeled  as  medium  1,  the  soft  coating  medium  is  labeled  as  medium  2, 
and  the  matrix  is  labeled  as  medium  3.  The  isotropic  material  constants  of  the  inclusion 
core,  the  coating  layer,  and  the  matrix  are  pi,  Fi,  Vi;  and  3,  respectively. 

The  radius  of  the  core  is  a.  The  lattice  constants  along  the  Yj  and  Y3  directions  are 
denoted  as  dj  and  1^3.  The  semimajor  and  semiminor  axes  of  the  ellipse  are  denoted  as 
^j^and  respectively.  In  the  unit  cell  (te,  i),  position  of  the  center  point  is  given  in  the 
global  coordinate  X2  =  and  ^3  =  ^3^.  For  convenience,  we  define  a  local  polar 

coordinate  system  (r,  8)  as  well  as  a  local  Cartesian  coordinates  (jr3,jC2)in  the  unit  cell 
with  jcs  =  rcosd  and  =  rsinS,  as  shown  in  Fig.  1(b).  Strong  anisotropic  properties  of 
the  metamaterial  along  the  X2  and  X^  directions  can  be  achieved  by  adjusting  dimensions 
of  the  major  and  minor  axes  of  the  coating  ellipse. 
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Fig.  6.1  (a)  An  anisotropic  elastic  metamaterial  made  of  cylinders  coated  with 
elliptical  soft  layer  in  a  matrix;  (b)  TheCfc,  0  element  (RVE). 


Due  to  geometric  complexity  of  the  microstructure,  analytical-based  methods  cannot 
be  applied  directly  for  the  determination  of  the  effective  dynamic  mass  density  [28].  In 
this  study,  a  numerical-based  method  (FEM)  will  be  adopted  [8].  Consider  the  unit  cell  in 
Fig.6.1(b),  time-harmonic  displacements  with  ce  =  2,3  are  prescribed  on  the 

exterior  boundary,  and  a  harmonic  analysis  of  the  unit  cell  is  performed  by  using  the 
finite  element  method,  from  which  resultant  forces  =  ^e^^on  the  exterior  boundary 
can  be  obtained.  The  effective  mass  density  tensor  can  then  be  numerically  defined  and 
determined  by 


where  =  OjOg  is  the  area  of  the  unit  cell.  To  demonstrate  that  the  mass  density  follows 
the  coordinate  transformation  law,  effective  mass  density  tensor  in  an  arbitrary  coordinate 
system  as  shown  in  Fig.2,  is  similarly  determined.  In  the  figure,  is  a 

rectangular  Cartesian  coordinate  system  in  principal  directions  and  (x  yj  is  an  arbitrary 
Cartesian  coordinate  system.  It  should  be  mentioned  that  the  unit  cell  in  Fig.  2  is 
specifically  suggested  for  simplification  of  numerical  calculation,  in  which  the  outer 
matrix  seems  to  follow  the  unit  cell  transformation.  This  is  because  the  outer  host  matrix 
is  mainly  functioned  as  an  effective  mass  in  the  calculation  of  the  effective  mass  density 
of  the  unit  cell.  Therefore,  determination  of  the  effective  mass  density  is  independent  on 
the  outer  matrix  geometry.  Based  on  the  numerical  results  of  the  unit  cell,  it  can  be 
proved  that  the  effective  mass  density  follows  the  coordinate  transformation  law  as 


3 


(6.2) 
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where  C  =  aas&,  S  =  sini,  and  £  is  the  angle  of  the  transformation.  Therefore,  the  effective 
mass  density  admits  the  second-order  tensorial  property. 


Table  6.1  The  microstructure  geometrical  parameters. 


Lattice  parameters 

Elliptical  coat 

Circle  coat 

=  1^3 

2Qmm 

20mm 

a 

3mm 

S.155mm 

7mm 

7.14mm 

^2 

5mm 

7.14mm 

Table  6.2  The  constituent  material  parameters. 


The  parameters 

Core:  lead 

Coating:  rubber 

Matrix:  epoxy 

Mass  density 

pi:  11,31  X  10^ 

p2:0,92  X  10^  kg/m^ 

pg:!,!!  X  10^  kg/m^ 

Young’s  modulus 

Ei:l,3 

E2;1,Ex  la^N/m^ 

E3:2,3EX  10’ AV 

Poisson’s  ratio 

Vi:  0,435 

I'j:  0,499 

V3:0,38 

Area 

A2:n(bib2  - 

.^3: 0^2 

Fig.  6.1  shows  the  effective  mass  density  of  the  elastic  metamaterial  along  a  direction 
£  =  30^  by  using  the  FEM  model,  where  flj-w  =  F  F  is  the  average 

static  mass  density  for  the  composite.  The  microstmcture  geometry  parameters  and  the 
constituent  material  parameters  are  given  in  Tables  6.1  and  6.2,  respectively.  From  the 
figure,  we  can  clearly  see  that  the  effective  mass  density  is  frequency-dependent, 

■s-ff  eff  eff 

anisotropic  ^  becomes  negative  around  the  resonant 

frequency.  For  the  metamaterial  with  the  current  microstructure,  two  different  resonant 
frequencies  along  and  Xg  directions  can  be  found  as  =  S62Hz  and  =  950Hz, 
respectively,  which  indicate  anisotropy  of  the  metamaterial.  The  significant  difference 
between  and  c?Ln  be  observed  when  the  frequency  is  close  to  the  resonant 
frequency.  It  is  also  interesting  to  note  that  the  mass  density  become  isotropic 
=  o)  when  the  frequency  is  much  less  or  much  larger  than  the  resonant 

frequency. 
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Fig.  6.1.  The  anisotropic  effective  mass  density  of  the  elastic  matatmaterial  predicted  by 
the  FEM  model  along  a  direction  S  =  30^ 

6.2  Fabrication  and  Testing  of  Anisotropic  Elastic  Metamaterial 

The  proposed  anisotropic  metamaterial  plates  were  manufactured  with  5=0  and 
5=90,  respectively.  The  elliptic  holes  were  first  cut  in  an  aluminum  (Al)  plate  (Type 
6061)  with  600mm  in  length,  33mm  in  width  and  3.175mm  in  thickness  using  a  computer 
numerical  control  (CNC)  machine.  The  major  and  minor  axes  of  the  elliptic  hole  are 
9.9mm  and  7.1mm,  respectively.  The  periodic  holes  (3X10)  are  arranged  in  a  rectangular 
lattice  with  lattice  constant  11  mm.  Then,  pure  lead  rods  (0=3. 175mm,  ESPI  Metals)  and 
epoxy  (p  =  1033k:g/m^,E  =  0.595GPa,u  =  0.37,  Crystal  Clear®  202,  Smooth-On  Inc.) 
were  filled  in  the  elliptic  holes  to  form  the  anisotropic  metamaterial  plate,  as  shown  in 


Eig.  6.2  A  fabricated  anisotropic  metamaterial  plate 
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Fig.  6.2.  The  center  positions  of  lead  rods  are  controlled  through  the  help  of  the 
additional  plate  with  exact  location  holes.  To  increase  anisotropy  of  the  metamtaerial, 
four  small  circular  holes  are  drilled  in  the  epoxy. 


6.2.1  Anisotropic  band-gap  validation 

To  investigate  the  band-gap  behavior  of  the  metamaterial,  the  transmission  property 
of  the  longitudinal  wave  propagating  along  the  X2  direction  is  first  studied 
experimentally.  The  set-up  of  the  transmission  experiment  is  shown  in  Fig.  6.3.  Two 
rectangular  piezoelectric  disks,  functioned  as  a  pair  of  actuators,  are  surface-bonded  to 
the  up  and  down  surfaces  of  the  plate  in  order  to  generate  longitudinal  wave.  Two 
circular  piezoelectric  disks,  functioned  as  a  pair  of  sensors,  are  used  to  receive  the 
transmitted  signals  after  the  metamaterial.  By  using  two  symmetric  sensors  on  the  plate, 
the  out-of-plane  components  can  be  eliminated  and  the  ratio  of  signal-to-noise  can  be 
significantly  increased  for  the  in-plane  component.  Fig.  6.4a  shows  the  experimental 
measurement  of  the  transmission  properties  of  the  metamaterial  at  5=0.  For  validation, 
the  exact  numerical  simulation  results  are  also  plotted  in  the  figure.  From  the  figure,  the 
band- gap  of  the  metamaterial  can  be  found  in  the  region  of  frequency  (17.2kHz, 
21.6kHz).  It  should  be  mentioned  that  the  ‘small  gap’  around  12.8kHz  from  the 
numerical  simulation  is  because  the  sensor’s  position  is  located  at  the  node  of  the 
standing  wave,  which  is  caused  by  the  band-gap  behavior.  Fig.  6.4b  shows  the 
experimental  measurement  of  the  transmission  properties  of  the  metamaterial  at  5=90..  As 
expected,  the  band-gap  of  the  metamateial  at  5=90.  is  different  from  the  metamaterial  at 
5=0  the  band-gap  frequency  is  in  the  frequency  regime  (23.8kHz,  31.6kHz).  However,  the 
experimental  results  can  accurately  capture  the  band-gap  behavior  of  the  metatmaterial. 


Fig.  6.3  Experimental  set-up  of  transmission  testing  of  the  metamaterial 
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(a)  (b) 


Fig.  6.4  Experimental  transmission  measurements  of  longitudinal  wave  in  the  anisotropic 
metamaterial  plates  for  6=0  (a)  and  5=90  (b) 


6.2.2  Determination  of  group  and  phase  velocities  of  the  metamaterial 

Analytically,  the  group  velocity  and  phase  velocity  of  the  longitudinal  wave 
propagating  along  the  X2  direction  in  the  metamaterial  plate  with  8=0  can  be  easily 
obtained  from  the  dispersion  relation.  Based  on  c,  =  —  ,  ,  the  group  velocity  can  be 
obtained  in  a  function  of  frequency  by  applying  central  difference  to  lowest  two 
longitudinal  branches  of  the  dispersion  relation  as  shown  in  Fig.  6.5(a).  The  phase 
velocity  can  be  obtained  from  as  shown  in  Fig.  6.5  (b). 


(a)  (b) 

Fig.  6.5  Group  velocity  (a)  and  phase  velocity  (b)  of  the  metamaterial  with  6=0 
obtained  from  the  lowest  two  longitudinal  branches  of  the  dispersion  relation. 
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In  order  to  experimentally  measure  the  group  and  phase  velocity  of  longitudinal 
wave  propagation  in  the  anisotropic  metamaterial  plate,  a  transient  experimental  set  up  is 
proposed  in  Fig.  6.6,  which  is  similar  to  the  set-up  of  the  harmonic  testing.  The  two 
rectangular  piezoelectric  patches,  functioned  as  actuators,  are  used  to  generate  in-plane 
longitudinal  wave.  Unlike  the  transmission  study,  two  pairs  of  symmetric  bonded 
piezoelectric  sensors  are  used  to  collect  the  in-plane  wave  signals.  Sensor  pair  1  is  before 
the  metamaterial  area  and  the  sensor  pair  2  is  after  the  metamaterial  area  as  shown  in  Fig. 
6.6.  In  order  to  capture  the  complete  wave  propagation  in  the  anisotropic  metamaterial 
plate  without  the  interference  of  the  reflected  waves,  the  total  length  of  the  anisotropic 
metamaterial  plate  is  extended  to  2.5m  by  welding  two  A1  plates  with  the  same  width  to 
the  two  ends  of  the  metamaterial  plate. 


Fig.  6.6  Experimental  set  up  for  group  velocity  and  phase  velocity  measurement  of  the 
metamaterial 

A  broadband  wave  signal  VCt)  =  ;4ill  -  cos (2^/01 5{n(2ff/0  is  chosen,  where  A  is  the 
amplitude  and  f  is  the  frequency.  Input  signals  with  f=13.5kHz,  15.5kHz  and  24.5kHz  are 
used  in  the  experiment  in  the  order  to  cover  the  interested  frequency  range  in  Fig.  6.5. 
The  continuous  wavelet  transform  (CWT)  based  on  the  Gabor  wavelets,  which  has  been 
demonstrated  to  be  very  useful  as  a  time-frequency  analysis  tool  of  wave  signals  in 
structural  health  monitoring  (SHM),  is  used  to  obtain  the  group  and  phase  velocity 
dispersion  curves.  The  continuous  wavelet  transform  (CWT)  of  a  given  signal  s(t)  is 
given  by: 


WTia,  b)  =  sCt^C— 5  Qt  (6.3) 

where  T)tr(t)  is  called  the  mother  wavelet  function  and  the  overline  indicates  the  complex 
conjugate,  a  and  &  are  known  as  the  scale  and  translation  parameters,  respectively.  The 
reciprocal  of  a  is  associated  with  frequency  and  &  is  related  to  time.  Gabor  function  is 
chosen  as  the  mother  wavelet  function  in  our  analysis  since  it  can  provide  a  better 
resolution  both  in  the  time  and  frequency  domain  than  any  other  wavelets  based  on 
Heisenberg  Uncertainty. 
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Fig.  6.7  shows  the  magnitudes  of  WT  coefficients  of  the  measured  output  signals 
collected  from  sensor  pair  1  and  sensor  pair  2  at  a  given  scale  o,  which  is  related  with 
certain  frequency  /.  The  arrival  times  for  the  output  signals  in  two  sensor  pairs  are 
determined  by  the  peak  locations  in  Fig.  6.7  (a)  and  (b),  respectively.  The  time 
difference,  At,  is  then  calculated.  Since  the  distance  between  sensor  pair  1  and  sensor  pair 
2  is  s’  =  215mm  and  the  length  of  the  metamaterial  area  is  =  110mm,  the  group  velocity 
as  a  function  of  frequency  can  then  easily  determined  as: 


^0 

- i- 

^3 


(6.4) 


where  is  the  group  velocity  of  the  matrix  material,  =  5400m/£  is  used  for  the 
aluminum  plate  at  low  frequencies. 

Figure  6.8  shows  the  experimental  measured  group  velocity  of  the  metamaterial  with 
5=0.  For  the  comparison,  analytical  solution  obtained  from  the  dispersion  relation  of  the 
metamaterial  in  Fig.  6.5  is  added  for  the  reference.  Good  agreements  between  the 
theoretical  result  and  the  experimental  measurement  can  be  found  in  the  most  interested 
frequency  regime. 
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_ Frequency  (KHz) _ 

Fig.  6.8  Comparison  of  the  group  velocity  from  the  numerical  simulation  and  experimental  testing. 


In  order  to  calculate  the  phase  velocity,  the  phase  angle,  of  the  WT 

coefficient  needs  to  be  considered.  At  a  given  scale  a,  the  corresponding  phase  angle  for 
the  output  signal  at  arrival  time  in  each  sensor  pair  is  determined  and  the  phase  difference 
between  the  received  signals  in  sensor  pair  1  and  sensor  pair  2  can  be  calculated. 
Then,  the  phase  velocity  is  obtained  in  a  function  of  frequency  using  the  following 
relation: 


(itC/)  F 


2wf 


rM  J 


where  At  is  the  arrival  time  difference  and  is  the  phase  velocity  of  the  matrix  material. 
Fig.  6.9  shows  the  comparison  of  the  phase  velocity  from  the  experimental  testing  and  the 
theoretical  result  from  the  dispersion  relation  of  the  metamaterial.  Similarly,  good 
agreements  can  be  observed  to  show  the  feasibility  of  the  proposed  experimental 
approach  for  determination  of  the  group  and  phase  velocities  of  the  elastic  metamaterial. 

Finally,  the  effective  dynamic  mass  density  of  the  anisotropic  metamaterial  can  be 
estimated  by  the  experimentally  obtained  phase  velocity  and  the  effective  stiffness 
matrix.  The  effective  stiffness  of  the  metamaterial  can  be  directly  determined  through  the 
unit  cell  analysis  and  is  found  that  it  is  not  dependent  on  the  wave  frequency.  Using  the 
equation  the  effective  dynamic  mass  density  of  the  anisotropic 

metamaterial  along  x  direction  can  be  experimental  determined  as  shown  in  Fig.  6.10. 
The  very  good  agreement  with  the  numerical  solution  can  be  also  observed. 
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Fig.  6.9  Comparison  of  the  phase  velocity  from  the  numerical  simulation  and 
experimental  testing 


Fig.  6.10  Comparison  of  the  mass  density  along  x2  direction  of  the  metamaterial  from  the 
numerical  simulation  and  experimental  testing 
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7.  DESIGN  OF  BROADBAND  ELASTIC  METAMATERIALS 


7.1  Introduction 

Although  many  metamaterials  have  been  proposed,  narrow  working  frequencies  of 
metamaterials  have  limited  their  applications  (Li  et  al.,  2009;  Capolino,  2009;  Cui  et  al., 
2010).  In  view  of  this,  there  are  several  papers  which  proposed  to  solve  the  narrow  band 
gap  problem.  However,  most  of  them  focused  on  optical  metamaterials;  there  are  just  a 
few  literatures  developing  acoustic/elastic  broadband  metamaterials.  Zigoneanu  et  al. 
(2011)  used  strong  anisotropic  acoustic  metamaterials  to  broaden  band  gap  up  to  500- 
3000  Hz.  Meng  et  al.  (2012)  applied  genetic  algorithm  to  achieve  the  optimized 
broadband  sound  absorption  800-2500  Hz.  In  this  study,  we  focus  on  the  design  of 
microstructures  to  improve  the  bandwidth  up  to  3370  Hz  starting  with  a  band-gap  lower 
bound  of  292  Hz.  Three  metamaterial  models  with  single  negativity  are  proposed,  and 
parametric  studies  as  well  as  comparisons  are  made  to  obtain  metamaterials  with  large 
band  gaps.  Moreover,  the  start  frequency  of  the  band  gap  is  in  hundreds  of  Hertz. 

7.2  Three  Models  for  Single  Negative  (SN)  Metamaterials 

Three  microstructures  with  dimensions  are  designed  for  metamaterials  as  shown  in 
Figs.  7.1  -  7.3.  The  material  properties  of  these  three  models  are  listed  in  Tables  7.1  - 
7.3,  respectively.  The  three  types  of  metamaterials  were  fabricated  with  a  3D  printer 
using  a  polymer  ink.  Figs  7.4  and  7.5  show  a  strip  of  the  Type  2  and  Type  3 
metamaterials,  respectively.  The  thickness  of  these  metamaterials  varies  over  the 
microstructure  in  order  to  achieve  the  equivalent  material  properties  listed  in  the  tables. 
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Unit :  mm 

Fig.  7.1  Type  1  SN  metamaterial  model 
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TABLE  7.1  Material 
Marked  area 

Blue 

Purple 

Green 

Orange 

Red 
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Unit :  mm 

Fig.  7.2  Type  2  SN  metamaterial  model 


parameters  for  the  Type  1  SN  microstructure 

Mass  density  p 

kg  / 

Young’s  modulus  E 

GPa 

Poisson’s  ratio  v 

1300 

5.0000 

0.49 

4000 

15.3846 

0.49 

1615 

6.2115 

0.49 

4000 

15.3846 

0.49 

1690 

6.5000 

0.49 

TABLE  7.2  Material  parameters  for  the  Type  2  SN  microstructure 


Marked  area 

Mass  density  p 

kg  / 

Young’s  modulus  E 

GPa 

Poisson’s  ratio  v 

Blue 

1300 

5.0000 

0.49 

Purple 

2000 

7.6923 

0.49 

Green 

1500 

5.7692 

0.49 
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, _ 40 _ , 

Unit :  mm 

Fig.  7.3  Type  3  SN  metamaterial  model 
TABLE  7.3  Material  parameters  for  the  Type  3  SN  microstructure 


Marked  area 

Mass  density  p 

kg  / 

Young’s  modulus  E 

GPa 

Poisson’s  ratio  v 

Blue 

1300 

5.0000 

0.49 

Purple 

2000 

7.6923 

0.49 

Green 

1500 

5.7692 

0.49 

Orange 

5200 

20.0000 

0.49 

Red 

5200 

20.0000 

0.49 

Fig.  7.4  Type  2  SN  metamaterial  manufactured  by  a  3D  printer 
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Fig.  7.5  Type  3  SN  metamaterial  manufactured  by  a  3D  printer 


In  the  Type  1  SN  metamaterial  model  shown  in  Fig.  7.1.,  the  two  green  horizontal 
beams  with  length  of  30  mm  and  width  of  1  mm  serve  as  two  vertical  springs  with  the 
desired  frequency  to  form  negative  Young’s  modulus  (NYM)  resonator.  Because  the 
rhombus  frame  is  only  used  for  transmitting  force,  the  width  of  frame  is  designed  as  1 .5 
mm  to  avoid  the  eigenmode  coupling  between  horizontal  beam  bending  and  frame 
stretching.  In  addition,  a  red  vertical  beam  serving  as  two  horizontal  springs  and 
supporting  triangular-shapes  masses  is  used  as  negative  mass  density  (NMD)  resonator. 
The  material  properties  are  chosen  to  broaden  the  bandwidth  by  combining  bandgaps  of 
NYM  and  NMD. 

Instead  of  using  beam  bending.  Type  2  SN  metamaterial  uses  frame  stretching  as  the 
stiffness  of  the  resonator.  As  shown  in  Fig.  7.2,  the  rhombus  frames  serve  as  two  vertical 
springs  supporting  purple  masses  to  achieve  NYM  local  resonance.  It  is  noted  that  Type  2 
SN  model  is  only  with  NYM.  The  geometry  and  material  parameters  are  chosen  from  a 
parametric  study  and  listed  in  Table  II. 

Type  3  SN  metamaterial  is  formed  by  adding  two  resonators  to  Type  2  SN  model.  As 
shown  in  Fig.  7.3,  the  two  red  vertical  beams  serve  as  two  sets  of  horizontal  springs  and 
support  orange  masses  to  achieve  NMD  local  resonance  when  they  move  in  phase.  In 
addition,  the  frame  stretching  resonator  and  the  symmetrical  mode  of  the  internal  masses 

(when  the  two  red  vertical  beams  move  180°  out-of-phase)  form  two  sets  of  NYM 
resonators  with  desired  frequencies.  The  band-gap  region  is  broadened  by  combining 
multiple  NMD  and  NYM  bandgaps. 

7.3  Mode  Shapes  and  Dispersion  Curves 

Considering  one-dimensional  longitudinal  harmonic  wave  propagation,  the 
displacement  can  be  expressed  as 

u(x,t)  =  (7.1) 

where  c  denotes  phase  velocity  and  q  wave  number.  By  selecting  two  locations  as  nodes 
(zero  displacement),  Eq.  (7.1)  can  be  treated  as  a  free  vibration  problem  of  a  finite  body 
with  the  fixed  end  condition.  With  this  concept,  we  constrain  the  horizontal  movements 
of  the  left  and  right  ends  of  a  100-unit-cell  long  strip  (so  the  two  ends  become  nodal 
points)  for  modeling  the  infinite  number  of  unit  cells  metamaterial  with  100  unit  cell 
periodicity.  The  result  of  the  free  vibration  analysis  yields  the  natural  frequencies  and  the 
corresponding  mode  shapes,  from  which  wave  numbers  are  obtained. 
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The  ABAQUS  commercial  software  is  used  for  the  finite  element  analysis.  Since  the 
beam  bending  serves  as  the  stiffness  of  resonator,  accuracy  for  beam  bending  estimation 
is  crucial.  For  element  type  selection,  the  incompatible  mode  element  (CPS4I)  is  applied 
due  to  its  efficiency  for  bending  calculation.  Four  CPS4I  elements  are  placed  along  the 
thickness  of  the  beam  to  achieve  element-number-convergence. 

Figs  7. 6-7. 8  show  the  dispersion  curves  for  Types  1-3  SN  metamaterials.  Each 
passing  band  corresponds  to  a  distinct  RVE  vibration  mode,  as  listed  in  Table  7.4  and 
Table  7.5. 

Dispersion  curves  for  the  Type  1  SN  microstructure  are  shown  in  Fig.  7.6,  the  NMD 
and  NYM  resonators  are  optimized  by  minimizing  the  passing  band  between  two 
adjacent  bandgaps  (the  narrow  passing  band  is  509.96—510.17  Hz).  It  is  noted  that  the 
other  narrow  passing  band  with  frequencies  463.35  —  463.47  Hz  in  Fig.  7.6  is  caused  by 
the  rotational  mode  of  the  horizontal  beams.  If  the  two  narrow  passing  bands  are 
neglected,  the  overall  band  gap  region  is  403-554  Hz,  firstly  caused  by  NMD  then  by 
NYM.  The  two  local  resonances  in  Type  1  SN  metamaterial  are  both  produced  by  beam 
bending. 


FIG.  7.6  Dispersion  curves  for  the  Type  1  SN  model 

Fig.  7.7  shows  the  dispersion  curves  for  Type  2  SN  metamaterial.  If  the  narrow 
passing  band  from  the  rotational  mode  ( 505.41-506.90  Hz)  is  ignored,  the  bandgap 
region  is  345  — 3097  Hz.  This  bandgap  is  caused  by  frame  stretching  which  produces  a 
much  broader  bandgap  than  the  beam  bending  mode. 
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FIG.  7.7  Dispersion  Curves  for  the  Type  2  SN  model 


TABLE  7.4  Mode  shape  of  the  Type  2  SN  model 


Frame  stretching  mode 


Rotational  mode 


340.85  Hz 


\rv 

501.93  Hz 


The  dispersion  curves  for  Type  3  SN  metamaterial  are  shown  in  Fig.  7.8.  There  are 
two  narrow  passing  bands  within  the  total  bandgap  region  of  292 -3694  Hz.  They  are  a 
rotational  mode  with  frequencies  752.40-755.24  Hz,  and  a  double  negativity  region  with 
frequencies  1064.8—1097.3  Hz.  The  overall  bandgap  region  is  caused  by  the  frame 
stretching  mode  (NYM),  the  symmetrical  mode  of  internal  masses  (NYM),  and  the  anti- 
symmetrical  mode  of  internal  masses  (NMD).  The  resonators  are  designed  in  order  to 
form  a  large  combining  band-gap  region.  It  is  noted  that  the  Type  3  SN  model  broadens 
the  bandgap  region  of  the  Type  2  SN  model  by  adding  two  extra  resonators. 
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FIG.  8.  Dispersion  curve  for  the  Type  3  SN  model 


TABLE  7.5  Mode  shape  of  the  Type  3  SN  model 


Frame  stretching  mode  Rotational  mode 


287.05  Hz 


746.23  Hz 


Internal  masses  symmetrical  mode  Internal  masses  anti-symmetrical  mode 


1062.0  Hz 


1089.8  Hz 


7.4  Representative  Volume  Element  Steady  State  Analysis 

In  addition  to  utilizing  dispersion  curves  to  obtain  the  bandgap  region,  we  also 
investigate  the  use  of  the  representative  volume  element  (RVE)  steady-state  analysis. 
Taking  Type  3  SN  model  for  example,  for  long  wave  motions,  the  frequency-dependent 
effective  modulus  of  the  metamaterial  can  be  obtained  by  applying  symmetrical  loading 
F^sinox  on  the  RVE  as  shown  in  Fig.  7.9(a).  The  displacements  at  side  edges  of  the 
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RVE  are  constrained  straight  and  with  the  corresponding  displacements  u  -  ±Mq  sin  ox . 
The  dynamic  stiffness  of  the  RVE  is  hence  given  by 


k 


=A 

2mo 


(7.2) 


It  is  noted  that  the  stiffness  of  the  unit  cell  defined  by  Eq.  (7.2)  is  proportional  to  the 
effective  Young’s  modulus  of  the  metamaterial. 

The  dynamic  stiffness  frequency  spectrum  is  shown  in  Fig.  7.9(b).  The  bandgap 
region  is  defined  as  the  frequency  range  within  which  k  has  negative  value,  i.e.,  294-866 
Hz  and  1064-3694  Hz. 

The  effective  mass  density  of  the  microstructures  can  be  obtained  from  the  dynamic 
response  of  the  RVE  loaded  as  in  Fig.  7.10(a).  The  displacements  at  side  edges  of  the 
RVE  are  set  to  be  equal  and  straight.  The  effective  mass  of  the  unit  cell  is  presented  by 


’^eff  =  - 


UgO) 


(7.3) 


The  effective  mass  frequency  spectrum  is  shown  in  Fig.  7.10(b).  The  bandgap  region 
is  defined  as  the  frequency  range  where  becomes  negative,  i.e.,  867-1099  Hz.  It  is 


noted  that  the  small  peak  in  Fig.  7.10(b)  at  frequency  755  Hz  is  caused  by  the  rotational 
mode,  which  can  also  be  observed  in  dispersion  curve  as  shown  in  Fig.  7.8. 

The  bandgaps  predicted  using  RVE  are  listed  in  Tables  VI  -  VII  together  with  the 
bandgaps  obtained  from  the  dispersion  curves  using  100  unit  cells  based  on  free  vibration 
analysis.  It  is  evident  that  the  RVE  steady-state  analysis  is  very  accurate  in  estimating  the 
bandgaps.  One  shall  also  notice  that  the  frequency  range  1065-1097  in  Table  VII  of  the 
dispersion  curve  approach  is  the  double  negativity  region  (band-gap  overlap  from  NMD 
and  NYM). 

The  reason  for  RVE  precisely  determining  band-gap  region  lies  in  the  accurate 
estimation  for  dispersion  curve  in  longest  and  shortest  waves.  When  the  wave  length  is 
much  larger  than  the  characteristic  length  of  RVE,  the  boundary  condition  for  each  unit 
cell  for  dispersion  curve  coincides  with  the  assumption  we  put  on  RVE:  the  side  edges  of 
RVE  have  equal  displacement  magnitude.  For  shortest  wave  length  in  each  branch  of 
dispersion  curve,  the  half  wave  length  is  exactly  the  same  as  the  characteristic  length  of 
one  unit  cell.  As  a  result,  each  unit  cell  in  dispersion  curve  has  the  equal  displacement 
magnitude  on  its  edge,  which  is  also  the  same  as  the  boundary  condition  we  use  for  RVE 
analysis.  It  shall  be  pointed  out  that  the  method  for  RVE  steady-state  simulation  has 
limitation  in  predicting  the  band-gap  region  for  vibration  mode  involving  non-straight 
side  edges  in  unit  cells. 


51 


Dynamic  stiffness  (N/m) 


{b) 


FIG.  7.9.  (a)  Testing  dynamic  stiffness  (effective  Young’s  modulus),  (b)  Dynamic 

stiffness  frequency  spectrum. 
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Fig.  7.10.  (a)  Testing  effective  mass  density,  (b)  Effective  mass  frequency  spectrum. 
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TABLE  7.6  Comparison  of  band  gap  region  prediction  for  Type  2  SN  model 
Method  Bandgap  Region  Estimation  (Hz) 

Dispersion  Curve  345  -  3097 

RVE  steady-state  analysis  347  -  3097 


TABLE  7.7.  Comparison  of  band  gap  region  prediction  for  Type  3  SN  model 

Method  Bandgap  Region  Estimation  (Hz) 

Dispersion  Curve  292  -  1065 

1097  -  3694 

RVE  steady-state  analysis  294  -  866  (From  RVE  NYM  test) 

867  -  1099  (From  RVE  NMD  test) 
1064  -  3694  (From  RVE  NYM  test) 


7.5  Design  Guide  for  SN  Metamaterials 

For  Type  1  SN  metamateraial,  the  wider  bandgap  in  high  frequency  range  is  easier  to 
produce  by  NYM  than  NMD.  Therefore,  a  smaller  mass  or  larger  stiffness  for  the  NYM 
resonator  is  used  to  broaden  the  bandgap  region  with  high  frequencies.  On  the  other  hand, 
for  the  NMD  resonator,  a  large  mass  and  large  stiffness  generally  result  in  a  broader 
bandgap,  and  the  resonator’s  stiffness  has  even  greater  influence  than  the  resonator’s 
mass.  However,  since  a  low-frequency  gap  is  desired  for  metamaterials,  large  resonator’s 
mass  is  still  recommended  for  the  NMD  model.  As  a  result,  the  appropriate  design  is  to 
use  large  mass  and  large  stiffness  for  the  NMD  resonator  to  forbid  low  frequency  waves, 
and  apply  small  mass  and  large  stiffness  of  the  NYM  resonator  to  block  high  frequency 
waves.  By  combining  NMD  and  NYM  bandgaps,  a  metamaterial  for  filtering  waves  with 
low  frequencies  in  wider  range  is  achieved. 

For  Type  2  SN  metamaterial,  a  large  mass  and  small  stiffness  of  the  NYM  resonator 
is  effective  to  shift  the  bandgap  region  to  low  frequencies.  However,  the  bandwidth 
shrinks  drastically  with  the  increment  of  the  mass  of  the  resonator.  As  a  result,  for  Type  2 
SN  metamaterial,  a  median  mass  of  the  NYM  resonator  is  recommended. 

For  Type  3  SN  metamaterial,  although  a  large  mass  of  two  NYM  resonators  (frame 
stretching  mode  and  internal  masses  symmetrical  motion  are  the  two  NYM  resonators) 
can  reduce  the  lower  bound  of  the  bandgap,  it  also  narrows  the  bandwidth.  Therefore,  a 
median  mass  of  the  NYM  resonator  is  recommended.  Trial  and  error  of  the  stiffness  of 
NYM  resonators  are  required  after  specifying  the  mass  of  NYM  resonator.  Since  the 
internal  masses  serve  as  NYM  and  NMD  resonators,  once  the  material  properties  of  the 
NYM  resonator  is  set,  the  NMD  resonator  is  determined. 
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7.6  Summary 

Three  types  of  SN  microstructures  are  designed  for  metamaterials  in  this  study. 
Parametric  studies  are  made  to  obtain  the  metamaterial  with  a  broad  bandgap  in  low 
frequency  range.  It  is  found  that  the  frame  stretching  mode  is  more  effective  than  using 
beam  bending  mode  in  designing  broadband  SN  metamaterials.  Moreover,  by  adding 
internal  mass  components,  the  bandgap  can  be  significantly  expanded. 

In  addition,  an  RVE  steady-state  analysis  to  determine  bandgap  regions  is  found  to  be 
quite  accurate  as  compared  with  that  based  on  dispersion  curves. 

The  broadband  metamaterials  proposed  in  this  study  are  fabricated  using  a  3D  printer. 

8.  REFERENCES 

Ao,  X.,  Chan,  C.T.,  2010.  Negative  group  velocity  from  resonances  in  two-dimensional 
phononic  crystals.  Wave  Random  Complex  20,  276-288. 

Ay  din,  K.,  Guven  K.,  2005.  Observation  of  negative  refraction  and  negative  phase 
velocity  in  left-handed  metamaterials.  Appl.  Phys.  Lett.  86,  124102. 

Boyd,  R.W.,  Gauthier,  D.J.,  2002.  Slow’  and  ‘Fast’  light.  Progress  in  Optics  43,  497-530. 
Boyd,  R.W.,  Gauthier,  D.J.,  2009.  Controlling  the  velocity  of  light  pulses.  Science  326, 
1074-1077. 

Brillouin,  L.,  1953.  Wave  Propagation  in  Periodic  Structures,  Dover,  New  York. 
Capolino,  F.,  2009.  Applications  of  Metamaterials.  CRC  Press. 

Cui,  T.J.,  Smith,  D.R.,  Liu,  R.,  2010.  Metamaterials:  Theory,  Design,  and  Applications, 
Springer. 

Chen,  H.,  Fung,  K.H.,  Ma,  H.,  Chan,  C.T.,  2008.  Polarization  gaps  and  negative  group 
velocity  in  chiral  phononic  crystals:  Layer  multiple  scattering  method.  Phys.  Rev.  B 
77,  224304. 

Cheng,  Y.,  Xu,  J.Y.,  Liu,  X.J.,  2008.  One-dimensional  structured  ultrasonic 
metamaterials  with  simultaneously  negative  dynamic  density  and  modulus. 
Physical  Review  B  77,  045134. 

Ding,  Y.,  Liu,  Z.,  Qiu,  C.,  Shi,  J.,  2007.  Metamaterial  with  simultaneously  negative  bulk 
modulus  and  mass  density.  Physical  Review  Letter  99,  093904. 

Drugan,  W.J.,  2007.  Elastic  composite  materials  having  a  negative  stiffness  phase  can  be 
stable.  Physical  Review  Letter  98,  055502. 

Fang,  N.,  Xi,  D.,  Xu,  J.,  Ambati,  M.,  Srituravanich,  W.,  Sun,  C.,  Zhang,  X.,  2006. 

Ultrasonic  metamaterials  with  negative  modulus.  Nature  Materials  5,  452-456. 
Feigenbaum,  E.,  Kaminski,  N.,  Orenstein,  M.,  2009.  Negative  dispersion:  a  backward 
wave  or  fast  light?  Nanoplasmonic  examples.  Optic  Express  17,  18934. 

Fok,  L.,  Zhang,  X.,  2011.  Negative  acoustic  index  metamaterial.  Phys.  Rev.  B,  83, 
214304 

Hu,  G.K.,  Huang,  G.L.,  Sun,  C.T.,  2011.  An  elastic  metamaterial  with  simultaneously 
negative  mass  density  and  bulk  modulus,  Appl.  Phys.  Lett.  98,  251907. 

Huang,  H.H.,  Sun,  C.T.,  2009.  Wave  attenuation  mechanism  in  an  acoustic  metamaterial 
with  negative  effective  mass  density.  New  Journal  of  Physics  11,  013003. 


55 


Huang,  H.H.,  Sun,  C.T.,  Huang,  G.L.,  2009.  On  the  negative  effective  mass  density  in 
acoustic  metamaterials.  International  Journal  of  Engineering  Science  47,  610-617. 

Huang,  H.H.,  Sun,  C.T.,  2011a.  Theoretical  investigation  of  the  behavior  of  an  acoustic 
metamaterial  with  extreme  Young’s  moduli.  J.  Mech.  Phys.  Solids  59,  2070-2081. 

Huang,  H.H.,  Sun,  C.T.,  2011b.  Locally  resonant  acoustic  metamaterials  with  2D 
anisotropic  effective  mass  density.  Philos.  Mag.  91,  981-996. 

Islam,  M.T.,  and  Newaz,  G.,  2012.  Metamaterial  with  mass-stem  array  in  acoustic  cavity. 
Appl.  Phys.  Lett.  100,  011904 

Jaglinski,  T.,  Kochmann,  D.,  Stone,  D.,  Lakes,  R.S.,  2007.  Composite  materials  with 
viscoelastic  stiffness  greater  than  diamond.  Science  315,  620-622. 

Lakes,  R.S.,  Lee,  T.,  Bersie,  A.,  Wang,  Y.C.,  2001.  Extreme  damping  in  composite 
materials  with  negative-stiffness  inclusions.  Nature  410,  565-567. 

Lakhtakia,  A.,  McCall,  M.W.,  Weiglhofer,  W.  S.,  2002.  Brief  overview  of  recent 
developments  on  negative  phase-velocity  mediums  (alias  left-handed  materials).  Int. 
J.  Electron.  Communications  56,  407-410. 

Lamb,  H.,  1904.  On  group-velocity.  Proc.  London  Math.  Soc.  1,  473-479. 

Lindell,  I.V.,  Tretyakov,  S.A.,  Nikoskinen,  K.I.,  Ilvonen,  S.,  2001.  BW  media  -  media 
with  negative  parameters,  capable  of  supporting  backward  waves.  Microw.  Opt. 
Technol.  Lett.  31,  129-133. 

Lee,  S.H.,  Park,  C.M.,  Seo,  Y.M.,  Wang,  Z.G.,  Kim,  C.K.,  2009a.  Acoustic  metamaterial 
with  negative  density.  Physics  Letters  A  373,  4464-4469. 

Lee,  S.H.,  Park,  C.M.,  Seo,  Y.M.,  Wang,  Z.G.,  Kim,  C.K.,  2009b.  Acoustic  metamaterial 
with  negative  modulus.  Journal  of  Physics:  Condensed  Matter  21,  175504. 

Lee,  S.H.,  Park,  C.M.,  Seo,  Y.M.,  Wang,  Z.G.,  Kim,  C.K.,  2010.  Composite  acoustic 
medium  with  simultaneously  negative  density  and  modulus.  Physical  Review 
Letter  104,  054301. 

Li,  J.,  Chan,  C.T.,  2004.  Double-negative  acoustic  metamaterial.  Physical  Review  E  70, 
055602. 

Li,  J.,  Fok,  L.,  Yin,  X.,  Bartal,  G.,  Zhang,  X.,  2009.  Experimental  demonstration  of  an 

acoustic  magnifying  hyperlens.  Nature  Materials  8,  931-934. 

Liu,  Z.,  Chan,  C.T.,  Sheng,  P.,  2005.  Analytic  model  of  phononic  crystals  with  local 
resonances.  Physical  Review  B  71,  014103. 

Liu,  Z.,  Zhang,  X.,  Mao,  Y.,  Zhu,  Y.Y.,  Yang,  Z.,  Chan,  C.T.,  Sheng,  P.,  2000.  Science, 
289,  1734  (2000). 

Meng,  H.,  Wen,  J.,  Zhao,  H.,  Wen,  X.,  2012.  Optimization  of  locally  resonant  acoustic 

metamaterials  on  underwater  sound  absorption  characteristics.  Journal  of  Sound  and 

Vibration  331,  4406-4416. 

Milton,  G.W.,  Willis,  J.R.,  2007.  On  modifications  of  Newton’s  second  law  and  linear 
continuum  elastodynamics.  Proceedings  of  the  Royal  Society  A  463,  855-880. 

Notomi,  M.,  2000.  Theory  of  light  propagation  in  strongly  modulated  photonic  crystals: 
Refractionlike  behavior  in  the  vicinity  of  the  photonic  band  gap.  Phys.  Rev.  B  62, 
10696. 

Schuster,  ]  A.,  1904.  An  Introduction  to  the  Theory  of  Optics.  Edward  Arnold,  London 

Shelby,  R.A.,  Smith,  D.R.,  Schultz,  S.,  2001.  Experimental  verification  of  a  negative 
index  of  refraction.  Science  292,  77-79. 


56 


Smith,  D.R.,  Padilla,  W.J.,  Vier,  D.C.,  Nemat-Nasser,  S.C.,  Schultz,  S.,  2000.  Composite 
medium  with  simultaneously  negative  permeability  and  permittivity.  Physical  Review 
Letter  84,  4184-4187. 

Pendry,  J.B.,  2000.  Negative  refraction  makes  a  perfect  lens.  Physical  Review  Letter  85, 
3966-3969. 

Veselago,  V.G.,  1968.  The  electrodynamics  of  substances  with  simultaneously  negative 
values  of  s  and  p.  Soviet  Physics  Uspekhi  10,  509-514. 

Wang,  L.  J.,  Kuzmich,  A.,  Dogariu,  A.,  2000.  Gain-assisted  superluminal  light 
propagation.  Nature  406,  277-279. 

Yao,  S.,  Zhou,  X.,  Hu,  G.,  2008.  Experimental  study  on  negative  mass  effective  mass  in  a 
ID  mass-spring  system.  New  Journal  of  Physics  10,  043020. 

Zigoneanu,  L.,  Popa,  B.I.,  Starr,  A.F.,  Cummer,  S.A.,  2011.  Design  and  measurements  of 
a  broadband  two-dimensional  acoustic  metamaterial  with  anisotropic  effective  mass 
density.  Journal  of  Applied  Physicsl09,  054906. 


9.  PUBLICATIONS 

•  H.  H.  Huang  and  C.  T.  Sun,  “Locally  Resonant  Acoustic  Metamaterials  with  2D 
Anisotropic  Effective  Mass  Density,”  Philosophical  Magazine,  Vol.  91,  No.6, 

2011,  pp.  981-996. 

•  R.  Zhu,  H.H.  Huang,  G.  L.  Huang,  and  C.  T.  Sun,  “Microstructure  Continuum 
Modeling  of  an  Elastic  Metamaterial,”  Int.  J.  Engineering  Science,  Vol.  49,  2011, 
pp.  1477-1485. 

•  H.  H.  Huang  and  C.  T.  Sun,  “A  study  of  Band-gap  Phenomena  of  Two  Locally 
Resonant  Acoustic  Metamaterials,”  J.  Nanoengineering  and  Nanosystems,  Nov. 
25,  2011,  pp.  83-91. 

•  X.  N.  Liu,  G.  K.  Hu,  C.T.  Sun,  and  G.  L.  Huang,  “Wave  Propagation 
Characterization  and  Design  of  Two-Dimensional  Elastic  Chiral  Metacomposite,” 
J.  of  Sound  and  Vibration,  330,  pp.  2536-2553,  2011 

•  X.N.  Liu,  G.  K.  Hu,  G.  L.  Huang,  and  C.T.  Sun,  “An  Elastic  Metamaterial  with 
Simultaneously  Negative  Mass  Density  and  Bulk  Modulus,”  Applied  Physics 
Letters,  98,251907,2011. 

•  H.H.  Huang  and  C.T.  Sun,  “Theoretical  Investigation  of  the  Behavior  of  an 
Acoustic  Metamaterial  with  Extreme  Young’s  Modulus,”  J.  Mechanics  and 
Physics  Solids,  doi:  10. 1016/j.jmps. 201 1.07.002,  2011. 

•  R.  Zhu,  G.  L.  Huang,  H.H.  Huang,  and  C.  T.  Sun,  “Experimental  and  Numerical 
Study  of  Guided  Wave  Propagation  in  a  Thin  Metamaterial  Plate,”  Physics  Letters 
A,  375,  ,2011,  pp.  2863-2867 

•  H.H.  Huang  and  C.T.  Sun,  “Continuum  Modeling  of  a  Composite  Material  with 
Internal  Resonators,”  Mechanics  of  Materials,  46,  2012,  pp.1-10. 

•  H.H.  Huang  and  C.T.  Sun,  “Anomalous  Wave  Propagation  in  a  One-dimensional 
Acoustic  Metamaterial  Having  Simultaneously  Negative  Mass  Density  and 
Young’s  Modulus,”  Journal  of  Acoustical  Society  of  America,  132(4), October 

2012,  pp.2887-95 


57 


•  R.  Zhu,  X.  N.  Liu,  Guoliang  Huang,  Hsin-Haou  Huang,  and  C.T.  Sun, 
“Microstructure  Design  and  Experimental  Validation  of  Elastic  Metamaterial 
Plates  with  Anisotropic  Mass  Density,”  Physical  Review  B,  Vol.  86,  No.  14, 
October  2012,  144307. 

•  K.T.  Tan,  H.H.  Huang  and  C.T.  Sun,  “Optimizing  the  Band  Gap  of  Effective 
Mass  Negativity  in  Acoustic  Metamaterials,”  Applied  Physics  Letters,  101: 
241902,  2012. 

•  Kota  Mikoshiba,  James  M.  Manimala  and  C.T.  Sun,  “Energy  Harvesting  Using  an 
Array  of  Multifunctional  Resonators,”  Journal  of  Intelligent  Material  Systems  and 
Structures,  Vol  24,  January  2013,  pp.  168  -  179. 

•  R.  Zhu,  X.  N.  Liu,  G.  K.  Hu,  ,  C.T.  Sun,  and  G.  L.  Huang,  “A  Chiral  Elastic 
Metamaterial  Beam  for  Broadband  Vibration  Suppression,”  J.  of  Sound  and 
Vibration,333,  2014,  pp.  2759-2773. 

•  K.T.  Tan,  H.H.  Huang  and  C.T.  Sun,  “Blast-Wave  Impact  Mitigation  Using 
Negative  Effective  Mass  Density  Concept  of  Elastic  Metamaterials,”  International 
Journal  of  Impact  Engineering,  Engineering,  Vol.  64,  2014,  pp.  20-29. 

•  James  M.  Manimala  and  C.  T.  Sun,  “Microstructural  Design  Studies  for  Locally 
Dissipative  Acoustic  Metamaterials,”  J  Applied  Physics,  115,  023518  (2014). 

•  James  M.  Manimala,  Hsin  Haou  Huang,  C.  T.  Sun,  Robert  Snyder  and  Scott 
Bland,  “Dynamic  Load  Mitigation  using  Negative  Effective  Mass  Structures,” 
Engineering  Structures,  Volume  80,  December  2014,  Pp.  458-468. 

•  Y.  Chen,  G.L.  Huang,  X.M.  Zhou,  G.K.  Hu,  and  C.T.  Sun,  “Analytical  coupled 
vibroacoustic  modeling  of  membrane-type  acoustic  metamaterials:  Membrane 
model,”  Journal  of  the  Acoustical  Society  of  American,  Vol.  136,  No.3, 
September  2014. 

•  Yongquan  Liu,  Xianyue  Su,  C.T.  Sun,  “Broadband  Elastic  Metamaterial  with 
Single  Negativity  by  Mimicking  Lattice  Systems,”  J  of  the  Mechanics  and 
Physics  of  Solids,  to  appear. 

•  R.  Zhu,  X.  N.  Liu,  G.  K.  Hu,  C.  T.  Sun  and  G.  L.  Huang,  “Negative  Refraction  of 
Elastic  Waves  at  the  Deep  Subwavelength  Scale  in  a  Single-Phase  Metamaterial,” 
Nature  Communication,  2014,  accepted. 

•  Y.  Chen,  G.L.  Huang,  and  C.T.  Sun,  “Band  Gap  Control  in  an  Active  Elastic 
Metamaterial  with  Negative  Capacitance  Piezoelectric  Shunting,”  ASME  Journal 
of  Vibration  and  Acoustics,  in  press. 


58 


